|Home | About | Journals | Submit | Contact Us | Français|
Short fragments of amyloidogenic proteins are widely used as model systems in studies of amyloid formation. Fragment 11–25 of the amyloid β protein involved in Alzheimer’s disease (Aβ11–25) was recently shown to form amyloid fibrils composed of anti-parallel β-sheets. Interestingly, fibrils grown under neutral and acidic conditions were seen to possess different registries of their inter-β-strand hydrogen bonds. In an effort to explain the microscopic origin of this pH dependence, we study Aβ11–25 fibrils using methods of theoretical modeling. Several structural models were built for fibrils at low and neutral pH levels and examined in short molecular dynamics simulations in explicit water. The models that display lowest free energy, as estimated using an implicit solvent model, were selected as representative of the true fibrillar structure. It is shown that the registry of these models agrees well with the experimental results. At the neutral pH, the main contribution to the free energy difference between the two registries comes from the electrostatic interactions. The charge group of the carboxy terminal makes a large contribution to these interactions and thus appears to play a critical role in determining the registry.
Amyloid fibrils are a class of large, filamentous aggregates that result from the self-assembly of a wide variety of peptides and proteins under appropriate thermodynamic conditions 1; 2. The recent vigorous scientific interest to amyloid fibrils is primarily due to their involvement in a host of progressive neurodegenerative diseases such as Alzheimer’s and Parkinson’s3. Additionally, amyloid fibrils also have well-defined and very specific functions in various living organisms4; 5 and they are rapidly becoming one of the most promising materials in nanotechnology6.
Understanding the causes of protein aggregation and characterizing the structure of fibrils are the first steps toward deciphering the role of these compounds in the context of the disease or designing their use in technological applications. As proteins related to disease are difficult to work with due to the high aggregation rate or associated safety risks, their fragments, which are easier to handle but still amyloidogenic, have emerged as a convenient model system for the general studies of amyloidosis7. Solid-state NMR experiments of Tycko et al8, have shown that a fragment containing residues 11 to 25 of the amyloid β peptide (Aβ11–25) involved in the Alzheimer disease forms amyloid fibrils. The fibrils display cross β-sheet structure9, which is the hallmark signature of amyloid, over the central hydrophobic cluster segment that comprises residues 17–21. Although experimental constraints were insufficient to solve the microscopic structure of the fibril, they provided strong evidence that the β-sheets are composed of β-strands in anti-parallel arrangement. Additionally, the registry of hydrogen bonds formed between β-strands was determined as 17 + k ↔ 20 −k at pH 7, where integers on both sides of the equation indicate residues in contact and k is an integer between −6 and 8. Remarkably, experiments repeated at an acidic pH 2 revealed that the β-sheets are still anti-parallel but their registry changed to 17 + k ↔ 22 − k, where k is an integer between −3 and 8. The shift induced by varying pH clearly indicates that the β-sheet registry in Aβ11–25 fibrils is determined by a sensitive balance among side-chain and main-chain interactions. However, with the atomic structure of the fibril absent, the exact mechanism by which the registry is encoded remains unknown.
Experimentally, solving the three-dimensional structure of an amyloid fibril presents a tremendous challenge10. In this paper, we follow an alternative route and try to determine the structure of Aβ11–25 fibrils using theoretical modeling and computer simulations. Computer simulations are widely and successfully used in structural studies of amyloid11; 12; 13; 14; 15; 16; 17. Our goal is to determine what interactions define the registry of the fibrils at low and neutral pH. First, we show through free energies calculations that β-strand dimers, which are the minimal building block of β-sheets, preferentially populate 17 + k ↔ 20 −k registry regardless of pH. At neutral pH, this indicates that the structure of the fibril is set at the lowest level of its hierarchical organization. For the acidic conditions, interactions between β-sheets making up a fibril appear to be important in determining the registry. Next, we build several two β-sheet models for Aβ11–25 fibrils under low and neutral pH and examine them in short, nanosecond-scale molecular dynamics simulations in explicit water. The models are scored according to their free energy estimated using generalized Born/solvent accessible surface area implicit solvation scheme. Conformations with the lowest free energy are selected as the true fibril structures for appropriate pH. It is seen that the registry of the generated structures at pH 2 and 7 agrees well with that determined experimentally. We find that under neutral conditions, a salt bridge between residues K16 and E22 plays a critical role in determining the structure of the fibril. As pH is lowered, the side chain of the glutamic acid becomes neutral thus reducing the impact of the salt bridge. Instead, hydrophobic packing and possible π – π interactions among aromatic residues emerge as a driving force in defining the structure of the fibril. Lastly, our models suggest that at pH 7, the major contribution to the free energy difference between registries 17 + k ↔ 22 − k and 17 + k ↔ 20 −k comes from the electrostatic interactions. By examining how these interactions are affected by changing the charge state of ionizable groups at varying pH, we were able to determine that the interactions involving charged C-terminus are responsible for a large part of the total free energy difference. This observation suggests that capping the C-terminus of Aβ11–25 with a neutralizing group may affect the registry of the fibril at pH 7, possibly reversing it to the 17 + k ↔ 22 −k configuration observed at pH 2. This hypothesis was tested and confirmed in computer simulations of the Aβ11–25 peptide with amidated carboxy terminus.
As pH changes from neutral to acidic, the charges undergoing transformation at ionizable sites in Aβ11–25 trigger a shift in the β-sheet registry. There are several residues in the primary sequence that change their charge state at acidic conditions. Negatively charged E11, E22 and D23 become neutral while neutral residues H13 and H14 turn positively charged. Additionally, carboxy group of the C-terminus becomes protonated and positively charged. At what level of the structural organization of the fibrils does the charge transformation affect the registry? Interactions among charged groups occur both locally, between nearest β-strands, as well as globally, between β-sheets or between protofilaments, as Figure 1(c) shows. Which of these interactions determine the registry? To address this question, we consider first the simplest possible model of a β-sheet, which is the one that consists of two peptides forming a β-strand dimer. An illustration of this model is shown in Figure 1(a)–(b). It serves as the elementary building block to construct larger β-sheets, allowing us to test the role of the local interactions. To estimate the relative occurrence rate of registries 17 + k ↔ 20 − k (R1) and 17 + k ↔ 22 ↔ k (R2) we computed their free energy difference. This is done in two steps. First, we use umbrella sampling18 combined with multiple histogram analysis method19; 20 to compute free energy profiles, or potentials of mean force (PMF), along the distance separating the geometrical centers of β-carbons located on two peptides that make up the dimer (see “Methods” for more details on the computations). The PMFs obtained at pH 7 and 2, G(r), shown in Figure 2, have a minimum at approximately rmin =5Å that corresponds to two β-strands being in contact and a plateau at large separations. The plateau appears at r = rp =8Å indicating the distance beyond which no additional cost is required to separate the peptides, or the distance beyond which they cease to interact. We can therefore interpret the difference in G(r) at the minimum and the plateau as the free energy cost of bringing two peptides together and forming a dimer, or the association energy. Figure 2 shows that this energy is approximately 30kCal/mol for R2 and 35kCal/mol for R1. Second, since the free energy of the state where peptides are completely dissociated (r > rp) is independent of the registry, we can use it as a common reference state to compute the relative free energy between two registries as the difference in their association costs. Our calculations show that the free energy of R1 is lower by approximately 5kCal/mol than that of R2, in agreement with the experimental observation. We conclude from these calculations that the registry of β-sheets in Aβ11–25 amyloid fibrils is set at the lowest level of their organization at neutral pH. Hydrogen bonds in anti-parallel β-strands together with local interactions among amino acid side chains, collectively favor R1 as the more stable configuration.
We repeated our calculations for acidic pH 2 and computed free energy profiles for R1 and R2. They look very similar to those shown in Figure 2 and therefore are not shown here. Surprisingly, we find that as at pH 7, the acidic conditions favor R1 over R2. The relative free energy actually increases and stands at about 10kCal/mol. This result shows that the local interactions among amino acid residues do not control the registry in Aβ11–25 fibrils.
Since the β-strand dimer can not explain why the registry changes at pH 2, it must be set at a different, higher level of fibril organization. As Figure 1(c) explains, the next structural element in fibrils hierarchy is β-sheets stacked against each other. It is not known how many β-sheets are contained in a protofilament of Aβ11–25. We focus on the minimal number, two, as the most convenient and computationally tractable model. It is common among fibrils formed by other proteins to find two β-sheets per protofilament21. Similarly to a single β-sheet, a pair of stacked β-sheets of arbitrary length can be constructed from an elementary building unit that consists of 4 peptides. A schematic representation of this unit is shown in Figure 3. As there is no experimental information regarding the mutual orientation of β-sheets in Aβ11–25 fibrils, we will consider all possible models. Figure 3 shows 4 different ways in which two anti-parallel β-sheets can be stacked. First, starting from two β-sheets perfectly superimposed, one β-sheet, for instance the left β-sheet, can be translated by the inter-sheet distance in the Y direction which is perpendicular to both the plane of the sheets and the fibril axis. We will call the resulting configuration C1. Three other configurations, C2, C3 and C4, can be constructed from C1 through three non-equivalent 180° rotations of its right β-sheet. One of these rotations is around the fibril axis while two others are around directions perpendicular to it. In the notations of Figure 3, these are rotations around X, Y and Z axes. Configurations C1 and C3 position the face of one β-sheet against the back of the other β-sheet and thus will be referred to as face-to-back. The other configurations, C2 and C4, maintain face-to-face contacts between the β-sheets. With respect to the mutual orientation of β-strands in the same β-sheet layer, C1 and C2 configurations can be classified as “parallel” and C3 and C4 as “antiparallel”. Therefore, an alternative and more intuitive classification for the stacking configurations is as follows: parallel face-to-back model for C1, parallel face-to-face for C2, antiparallel face-to-back for C3 and antiparallel face-to-face for C4. Our stacking configurations C1–C4 correspond to the symmetry classes 7 and 8 of the cross-β sheet classification introduced by Eisenberg et al22.
In addition to the registry, R1 or R2, and stacking, C1–C4, one more variable is needed to uniquely define the structure of a two-β-sheet protofilament. That variable describes the mutual shift of the two β-sheets along the lateral direction of the fibril. Using the notations of Figure 3, this direction is along X axis. The mutual shift of β-sheets determines which contacts are formed at the interface, or the sheet-to-sheet registration, and therefore is an important characteristic of a fibril model. To define the inter-sheet displacement, we will use notations Dk, where k is an integer between −7 and 7 that indicates that in the top layer of our four-peptide β-sheet model (see Figure 3), residue F19 of the left β-strand is in contact with residue 17+k of the right β-strand. This allows us to adopt Ri-Cj-Dk notation, where i=1,2, j=1–4, k=−7,7, that unambiguously defines all possible configurations of two stacked β-sheets of a particular registry. As an example, Figure 4 shows C3 configurations of registry R1 and inter-sheet displacement D1 and D−1, R1 – C3 – D1 and R1 – C3 – D−1. Both structures are in the antiparallel face-to-back orientation. The first includes inter-β-sheet contacts between F19 and V18 while the second between F19 and K16.
How β-sheets are stacked and displaced with respect to one another in Aβ11–25 fibrils is not known. We therefore attempt to determine stacking configuration C and displacement D theoretically. In what follows we consider all possible C–D combinations and select the most likely candidates for a fibril model based on the results of computer simulations. One important condition that we use in this procedure is that fibril models do not contain uncompensated charged residues buried in their interior. In large proteins and protein complexes, charged side chains are found either on the surface when exposed to solvent, or engaged in salt bridges when in the complex’s interior23. We determine the most stable model for normal and low pH separately.
We start our analysis with the configuration C3 in registry 1 at pH 7. As Figure 4 shows, it displays a tendency to bury uncompensated charged groups inside the fibril. In R1 –C3 –D1 for instance, charged lysines, glutamic and aspartic acid side chains are located deep inside the fibril disrupting the hydrophobic core. It is easy to see that any meaningful β-sheet displacement can not relieve this constraint. By shifting the right β-sheet up or down by two residues, it is possible to create a favorable contact between D23 and K16 in one layer, as in configuration R1 – C3 – D−1. This pair of residues, however, remains uncompensated in the other layer, as seen in Figure 4. Additionally, negatively charged side chains of E22 still remain inside the fibril. We conclude from this analysis, that configurations C3 are not suitable for a fibril model. Conformations in C1 and C2 stacking configurations are shown in Figure 5(a)–(b). They are seen to either (i) place uncompensated charges inside the fibril as in R1 – C1 – D3, or (ii) create contacts between charged residues of the same sign as in R1 – C2 – D2. Like in C3 configurations, these conflicts can not be relieved by any meaningful displacements of the β-sheets and, therefore, C1 and C2 configurations will not be considered further.
We find that the only configuration that can satisfy the constraints on charged residues is C4. Two conformations of this configuration are shown in Figure 5(c)–(d). In conformation R1 – C4 – D0, residues K16 are in contact with residues E22 in the bottom layer. The other charged residues in this layer, D23 and E11 are found on the surface of the model. In the top layer, residues D23 and E11 are in the interior of the fibril but sufficiently close to the edges such that they may get exposed to solvent. Residues K16 and E22 are found on the surface of the model. The stability of R1 – C4 – D0 model was tested in short nanosecond-scale simulations in explicit solvent. The initial conformations were constructed manually as discussed in more detail in the “Methods” section. We used 8 Aβ11–25 peptides per β-sheet, or 16 peptides per protofilament, to construct our model. It is known that fibril assembly follows the classical nucleation-growth scenario24. To be representative of an intact fibril, theoretical models need to consider number of peptides greater than the size of the critical nucleus. In studies of amyloid formation, it is not uncommon to observe relatively small nuclei, some of them containing 6 monomers or less25. Furthermore, peptide concentration, temperature and other parameters define critical nucleus and thus can be used to control its size. Although we do not know what the nucleus is for Aβ11–25 fibrils grown under the reported conditions8, we assume that it is possible to reduce its size to fewer than 16 monomers by selecting appropriate thermodynamic conditions.
Our protofilament model was solvated in explicit solvent, making sure that there were no water molecules in its interior. After minimization, the system was equilibrated for 20ps with all the peptide atoms fixed. The equilibration phase was followed by an 8ns production run where both peptides and water were fully unconstrained. The initial configuration used in our simulations is shown in Figure 6(a),(c). It shows β-sheet structure, according to DSSP classification26, extending over the maximum number of residues permitted under registry 1 and encompassing residues H13 to V24. This initial secondary structure survives 8ns of simulation almost intact with a slightly reduced β-sheet content of H13 and V24 toward the end of the trajectory but still populated more than 90% of the time. More detailed structural analysis of β-sheets observed in our simulation is presented in Fig. 1S (Supporting Information). An important feature that emerges in these structures is a pattern of alternating hydrogen bonds between C- and N-termini of neighboring β-strands that propagate throughout the length of the fibril.
At the end of the 8ns simulation, the initial fibril configuration undergoes a noticeable conformational change. Two projections of the final conformations are shown in Figure 6(b) and (d). Among other structural features, the model develops a twist, which is a common characteristic of amyloid fibrils27 and is believed to optimize the hydrogen bonds, side chain stacking and electrostatic interactions. Another apparent structural change is that the fibril begins to open up at its edges, as seen in Figure 6(d). This results in a less compact structure at its boundary with water, where the void created by diverging β-sheets is filled with solvent. As Figure 6(e) shows, depicting the microscopic details of the fibril edge, water penetrate deeply into the fibril, completely solvating the side chains of residues E11-H14 and D23-G25. This is a direct consequence of many of these residues being either hydrophilic or charged and therefore favorably interacting with water. A cluster of residues between K16 and E22 in the central part of the sequence are seen to be completely protected from water. As Figure 6(f) shows, this cluster is centered around a group of strong contacts between residues F19 and L17 and between V18 and F20, which form an extensive and insoluble hydrophobic complex in the middle of the fibril conformation. The presence of such hydrophobic clusters, or cores, is a hallmark feature of amyloid structure7.
The hydrophobic core is flanked on both sides by salt bridges between K16 and E22, which protect it from water. A snapshot of a representative water configuration is shown in Figure 6(e). The salt bridges appear to be forming a “wall” made of positive and negative charges that benefit from the presence of water. The hydrophobic residues behind the wall never come in contact with water and thus avoid unfavorable interactions. Additionally, the side chain of Q15, which is aligned with A21 on the same layer of β-sheet (see Figure 5(c)), is seen to form a hydrogen bond with the side chain of K16 located on the other layer. These interactions appear to make the wall even stronger.
To test how credible the results of our R1 – C4 – D0 simulations are we designed a control run in which a fibril model is studied that is expected to be unstable. We considered the R1 – C4 – D2 model, shown in Figure 5(d), which can be constructed from R1 – C4 – D0 by shifting one of the β-sheets up or down by two residues in the lateral direction, and which contains an internal conflict in the location of charged residues. It is seen that the constraints on charged residues are all satisfied but K16 and E22 are not aligned. Residues E11 and D23 are located near the edge of the fibril, and thus may get exposed to solvent. Additionally, the model displays what seems an uninterrupted hydrophobic core, possibly stabilized by favorable aromatic interactions of phenylalanine side chains. To test whether R1 – C4 – D2 conformation is stable on nanosecond time scale, molecular dynamics simulations were performed using the same protocol as for R1 – C4 – D0. We find that the initial β-sheet structure quickly collapses after only 2ns of the simulation time. A snapshot of the collapsed conformation is shown in Fig. 2S (Supporting information). The collapse appears to be prompted by the uncompensated charges of K16 side chains interacting with backbone oxygens of adjacent β-strands, resulting in the creation of a hole-like defect in the β-sheet wall. The defect is accessible to water molecules which were seen to flow freely into the interior of the fibril. Our control simulation thus indicates that R1 – C4 – D2 conformation is unstable, justifying the use of the short-time molecular dynamics simulations in this work.
To construct fibril models in registry 2 at pH 7 we followed the same procedures as described above for registry R1. Four different stacking configurations were considered for two β-sheets. Figure 7 shows the only two configurations, R2 – C4 – D4 and R2 – C3 – D3, that either expose charged groups to solvent or maintain them in salt-bridge pairs. Conformation R2 – C4 – D4shares with the R1 conformation R1– C4 –D0 its distinctive K16-E22 salt bridge. In contrast, conformation R2– C3 – D3 features K16-D23 salt bridges. Both conformations were tested in 8ns molecular dynamics simulations, where they produced stable trajectories maintaining the initial β-sheet structure.
To differentiate among multiple conformations that appear to be equally stable, we evaluate their relative free energy. For compounds as complex as fibrils, it is impossible to compute free energy differences exactly, as we did for the β-strand dimer. We therefore resort to approximations here. A combination of the Generalized Born (GB) model for the electrostatics part of the solvation free energy with the solvent accessible surface area (SA) model for the non-polar part are used to estimate free energy of fibril models. More detailed description of the GB/SA scheme is given in the “Methods” section.
Since we can not fully equilibrate the simulated systems globally, we require that local equilibration occurs before estimating the free energy. Provided that the initial configurations are sufficiently close to the equilibrium states, the local estimates obtained in simulations provide a good representation of the global free energy. To monitor how fast equilibration occurs, we examine total free energy as a function of time. Two examples of such time dependences are shown in Figure 8(a) for conformations R1 – C4 –D0 and R2 – C4 – D4. It is seen that the two curves share common behavior: they start at a relatively high value (reflecting that the initial state is removed from equilibrium) and rapidly fall off to a plateau where they stay throughout the remainder of the simulation. The length of this rapid equilibration phase depends on the fibril model and stands at ~1ns for R1 – C4 – D0 conformation and at ~2ns for R2 – C4 – D4 model. To determine the free energy of locally equilibrated states, G is averaged over the plateau area. As Figure 8(a) shows, the averaging time varies with the specific fibril model. A total of 100 configurations are selected uniformly over the averaging interval to estimate free energies and their errors. The resulting free energy traces, broken into components, are shown in Figure 8(b) for the R1– C4 – D0 model. It is seen that the electrostatics component completely dominates the total free energy, being nearly two orders of magnitude greater than the other two parts.
Free energy estimates for the R2 family of structures, R2 – C4 – D4 and R2 – C3 – D3 (see Figure 7) show that the latter is higher in free energy by approximately ΔG= 30kCal/mol. We will therefore focus on R2 – C4 – D4 structure as the most stable model of amyloid in registry R2. Figure 8(a) compares the free energy of this structure with the free energy of R1 – C4– D0 model, which is the most likely fibril model in registry R1. It is clear that registry R1 is more stable at pH 7 and that the free energy difference between the two competing structures is statistically significant. Broken into parts, ΔG = G2−G1= 143±16 kCal/mol is composed of ΔΔGel = 187 ± 8kCal/mol, ΔUvdw = 28±6kCal/mol and ΔΔGnp= −72±52kCal/mol. The statistical errors in these expressions were estimated assuming that the conformations over which the averaging is carried out are statistically independent. This obviously is not the case in our nanosecond-scale simulations. To get a better estimate of the errors, we carried out additional two simulations for R1– C4– D0and R2 – C4 – D4 structures, using slightly different initial structures. These simulations (data not shown) indicated that the total free energy difference ΔG is reproducible, with the estimated error increasing slightly from 16kCal/mol originally to 20kCal/mol in additional tests. We also found that the individual free energy components are subject to larger fluctuations than their sum, which is understandable.
We estimate that the error in the electrostatics part ΔΔGel increases to close to 30kCal/mol, from the value of 16kCal/mol suggested originally. But even with this increase, the free energy difference between R1 and R2, ΔΔG el = 187kCal/mol, remains statistically meaningful. Based on Figure 8(a), we conclude that the lowest free energy conformation for Aβ11–25 fibril at pH 7 contains β-sheets in registry 1.
In contrast to neutral pH, at the acidic pH 2 the carboxylic groups of glutamic acid, residues E11 and E22, aspartic acid, residue D23, and the C-terminus become neutral. At the same time, the histidine residues H13 and H14 become protonated and join K16 and N-terminus as the only charged groups in the peptide. Importantly, since these groups carry the same charge +1, they can not form salt bridges inside fibril’s interior as discussed in the previous sections for pH 7. Consequently, one of the constraints we use to construct fibril models is eliminated. The second constraint, however, requiring that charged residues be found in contact with solvent, remains intact.
We examined all four stacking families C1–C4 with registry 1 and found only 4 conformations that may potentially satisfy this constraint. These conformations, R1 – C4 – D0, R1 – C4 – D2, R1 – C1 – D5 and R1 – C1– D7, are depicted schematically in Figure 9. All four models were simulated for 8ns following the protocol described in the previous sections. Free energy analysis indicates that the most stable states are R1– C4– D2 and R1– C1 – D5 whose free energies are statistically indistinguishable. The free energy of two other states is higher by 44kCa/mol for R1– C4– D0 and by 56kCal/mol for R1– C1 – D7. That there are two states with the same free energy seems surprising at first, given the large apparent differences between these states, compare Figure 9(b) and (c). A closer examination, however, reveals that R1– C4– D2 and R1– C1– D5are related. These conformations can be transformed one into another by a shift of one of the two β-sheets (see Figure 3) by one layer along the fibril axis (the direction of the shift is unimportant). As Figure 9 shows, all that results from the shift is a subtle rearrangement of the interfacial contacts, which are mostly hydrophobic in nature. Note that such a shift would not be possible in R1 – C4 – D0 at pH 7, as it would entail the disruption of the K16-E22 salt bridges. Whether R1 – C4 – D2 and R1 – C1 – D5conformations truly have identical free energies or they only appear to do so because of our indiscriminate treatment of hydrophobic interactions modeled through solvent exposed surface area is not clear. But within the validity limits of our model, these conformations are identical. In what follows, we will limit our analysis to only one of them, R1 – C4 – D2.
We find that a total of 7 registry R2 models can satisfy the constraint on the charged residues. Two of these models, R2 – C4– D0 and R2– C4– D2are shown in Figure 10 while five others, R2 – C4 – D4, R2– C1– D5, R2– C1– D3, R2– C1 – D1and R2– C3– D3 are shown in Fig 3S (Supporting information). Compared to R2 – C4 – D0, conformation R2 – C4 – D2 seems to offer better hydrophobic contacts, possibly aided by π − π interactions among side chains of F19 and F20, but it may be more restrictive in terms of access of K16 side chains to water. All seven models produced stable 8ns trajectories in which they retained the initial β-sheet structure and displayed other characteristics of amyloid fibrils.
Free energy analysis of the simulated structures suggests that pairs of conformations R2 – C4– D0 -R2– C1– D5, R2 – C4 – D2 - R2– C1– D3 and R2– C4 – D4- R2 – C1 – D1have identical free energies, within statistical errors. These conformations are related by the same shift transformation as discussed for the registry 1 models. We will therefore limit our discussion to the C4 family only. We find that among conformations R2 – C4 – D0, R2 – C4 – D2 and R2 – C4 – D4, the latter is approximately 100kCal/mol higher in free energy than the former two, and will not be further analyzed. The free energies of R2– C4– D0 and R2– C4– D2 states appear to be indistinguishable within statistical errors. On the whole, this is not surprising, given that the two structures are closely related and that the main forces holding two β-sheets together in these models are non-specific hydrophobic interactions. However, a more detailed analysis of the final conformations generated for these structures, indicated that one of their β-sheets is rotated with respect to the other β-sheet by about 20°. A snapshot of the rotated conformation is shown in Fig. 4S (Supporting information). The rotation is not consistent with the fibril architecture. To test whether it could affect the structure of simulated models to any significant extent, or change the ranking of free energy, a model that contains 16 peptide chains per β-sheet (32 chains in total) was investigated. Initial conformations for R2 – C4– D0 and R2 – C4 – D2 models were built and simulated in exactly the same manner as discussed in the previous sections. The simulations were run for 6ns and were repeated twice to test the convergence. We find that the β-sheet rotation observed in the 16-peptide model is no longer present in the 32-chain simulation. But the states R2– C4 – D0 and R2– C4 – D2 still display essentially the same free energy. It appears, therefore, that either the Aβ11–25 fibrils are characterized by a genuine structural polymorphism at pH 2, or our free energy model is not accurate enough to detect the difference between two alternative states. We will assume the former is the case and focus our analysis on one of the equivalent structures, R2– C4– D0.
In addition to R2 – C4 – D0 model, a free energy estimate was generated for a 32-chain model of the most stable R1 structure, R1 – C4– D2. This allows for the direct comparison of the free energy between the two competing structures, which is shown in Figure 11. It is seen that the time dependences closely resemble those discussed in Figure 8: rapid fall-off followed by a plateau. Averaging was performed over the last 40 structures selected for the analysis, as shown in the inset of Figure 11. It is clear that registry 2 is the more stable of the two structures. The corresponding free energy difference ΔG= G2 − G1= 118± 27 kCal/mol can be decomposed into contributions of ΔΔGel= −146±15kCal/mol, ΔU vdw= 156±10kCal/mol, and ΔΔGnp= −129±3kCal/mol. An analysis of the free energy difference obtained using the 16-chain model yields similar results.
Snapshots of the last conformation of our R2– C4–D0 simulation are shown in Figure 12. The model has many structural properties, such as twist and hydrophobic core, in common with the model R1– C4– D0 that was derived for pH 7. The initial β-sheet structure of the two β-sheets survives the simulation almost intact. In the last 1ns, residues Q15-V24 were seen in the β-sheet state more than 90% of the time. The number of residues in the β-sheet configuration is by 2 residues less than at pH 7. This is fully consistent with the registry of the R1 – C4– D0 model. More detailed analysis of the β-sheet structure observed in the simulation is presented in Figure 5S (Supporting information). The map of hydrogen bond probabilities, shown in this figure, is consistent with the secondary structure analysis, indicating strong bonds involving residues Q15 to V24. Importantly, there is no hydrogen bond observed between C- and N-termini as was discussed for the R1 – C4 – D0 structure at pH 7.
Like model R1 – C4 – D0 at pH 7, the R2 – C4 – D0 model at pH 2 displays completely solvated edges. The extent of water penetration is illustrated in Figure 12(b). The only part in the fibril shielded from solvent contains tightly packed side chains of L17-E22. These residues form the hydrophobic core, which is highlighted in Figure 12(c). Interestingly, we find that the core is completely dry around F19-L17 and F20-E22 contacts, while Q15-A21 and V18-V24 contacts allow partial solvation. The map of inter-β-sheet contacts shown in Figure 5S confirms this conclusion. Fulfilling the design requirement, the charged side chains of K16 are fully solvated in water. A detail of the solvated configuration is shown in Figure 12(d). Water is seen to freely flow around charged amino groups of K16. In the space formed between K16 and V18 and between D23 and A21, a string of highly ordered water molecules is seen, forming a 1-D water wire. The significance of this observation is not clear at the moment.
The main structural difference between Aβ11–25 fibrils grown at different pH concerns the registry of their β-sheets. Specifically, at pH 7, residues 17+k of one β-strand are aligned with the residues 20-k of the next β-strand thus forming registry 17 + k ↔ 20 −k (R1). At pH 2, the registry changes to 17+ k ↔ 22 −k (R2), as shown in Figure 1(a)–(b). As the pH is lowered, several ionizable groups change their charge state thereby affecting interactions among β-strands. The resulting perturbation of the balance in the electrostatic interactions leads to the registry shift. It is possible that local interactions operating between two adjacent β-strands in a β-sheet trigger this structural transition. To test whether this scenario takes place, we evaluated relative free energy of forming a β-strand dimer in registry R1 and R2. We find that registry R1 is more stable at both neutral and low pH. As Figure 1 shows, R1 allows for two more inter-strand hydrogen bonds than registry R2, and thus is favored by hydrogen bond interactions. Our results suggest that β-sheet registry in Aβ11–25 fibrils is determined by long-range or global, rather than local, interactions among peptides. An important consequence of this assertion is that at acidic pH, small oligomeric species that precede fibril formation will be trapped in the wrong registry. This may have kinetic consequences for the fibril growth, as oligomeric intermediates would have to be dissociated first before being adsorbed onto fibrils, thus slowing down the overall process. Our simulations therefore suggest that, in addition to the structure of fibrils, varying pH also affects fibrillization pathways.
To try to explain how global interactions among peptides are involved in setting up the registry, we constructed amyloid models composed of two stacked β-sheets, corresponding to the next level in the structural organization of fibrils above β-strand dimers. Simulations of these models show that the most stable structure at pH 7 has the antiparallel face-to-face stacking configuration with the F19 residues on one β-sheet forming contacts with L17 residues on the other β-sheet. The registry of this structure is 17 k ↔ 20 − k (R1), which agrees well with the experiment. At pH 2, structures in the same stacking configuration and featuring contacts of F19 with either L17 or F20 emerge as most stable. Again, the predicted registry 17 + k ↔ 22−k (R2) is in good agreement with experiment, lending support to the employed computational method and the derived structural models. With the microscopic conformations of Aβ11–25 fibrils available (to the extent that our model is correct), one can begin to gain insights into what interactions control their assembly.
We find that at both low and normal pH, the derived structures exhibit large clusters of tightly packed non-polar residues, or a hydrophobic core, which is a distinguishing characteristic of amyloid fibrils28. While readily solvating fibrils’ edges, water is expelled from their cores in our simulations, in accordance with the hydrophobic effect. The resulting interface between β-sheets is dry in the central part of the fibril. This feature is in common with the steric zipper model of amyloid fibrils, suggested on the basis of crystallographic studies of short amyloidogenic sequences22. Unlike this model, however, we do not see an interdigitation of residues at the interface. Hydrophobic residues from different β-sheets are seen to form dehydrated contacts but we have not observed a pattern of penetration of the side chains on one β-sheet into the interspace between the side chains of the opposite β-sheet, consistent with the steric zipper. It is not clear from our results whether the Aβ11–25 fragment does not form steric zippers or our simulations were not sufficiently long to observe their formation.
At neutral pH, the charged residues E22 and K16 play a key role in determining the structure of the fibril. In registry R1, a configuration that allows for E22-K16 salt bridges is seen as the only stable structure that survives a nanosecond-scale simulation intact. For registry R2, two configurations exists, but the one with the E22-K16 contacts appears to be more stable. As the pH is lowered to acidic, the charged residue E22 becomes neutral thus removing the beneficial effect of the salt bridge. As a consequence, the role of the electrostatic forces diminishes while other types of interactions, such as hydrophobic and dispersion forces, gain influence. To quantify the role of different interactions we examine free energy difference ΔG = G2 − G1 between R2 and R1 structures. Three different components of ΔG are presented in Table 1 at pH 2 and 7. It is seen that at pH 7, the more stable R1 structure is favored by electrostatic forces as well as van der Waals interactions. The stabilizing effect of the electrostatic interactions agrees with the recent simulations of Thirumalai et al on an unrelated Aβ fragment29. Registry 2 conformation is seen to offer better hydrophobic interactions. At pH 2, contributions from different free energy components are, for the most part, consistent with those observed at pH7. Van der Waals forces favor the R1 structure while hydrophobic forces – R2 conformation, reflecting that these short-range interactions are mostly defined by the geometry of the particular registry. Electrostatic forces, on the other hand, change their preferred structure depending on the pH of the environment. This is not surprising since a change in pH affects the charge state of ionizable groups.
Although the electrostatic forces favor the most stable state at both neutral and low pH, their relative contribution differs in these two environments. At pH 7, ΔΔGel completely dominates the total free energy ΔG. It is more than 6 times greater than the van der Waals contribution, and more than 2 times greater than the non-polar term. To a large extent, ΔΔGel determines the registry at pH 7. In contrast, all three free energy components, ΔΔGel, Uvdw and ΔΔGnp have comparable magnitudes at pH 2. Consequently, the registry is determined by a delicate balance among various long-range interactions. The fact that ΔΔGel is so large at pH 7 makes it a suitable target for rational manipulations, designed to alter the total free energy difference and thus shift the balance between fibrils of different registries. A striking result from Table 1 is that ΔΔGel changes its sign when the pH is lowered from neutral to acidic. We know that the pH change alters the charge state of six ionizable groups in the primary sequence of Aβ11–25, C-terminus, E11, E22, H13, H14 and D23. We also know that the end result of this alteration is ΔΔGel that favors registry 2 instead of registry 1. Is it possible that some charge groups in the six involved residues contribute to ΔΔGel more than others? This would be a pertinent question to ask in the context of rational fibril design. Altering the charge state of key residues through mutations would allow one to design a sequence that forms registry 2 fibrils at pH 7 instead of registry 1 fibrils. Answering this question would also clarify the role of various residues involved in defining the registry in Aβ11–25 fibrils.
With the trajectories of fibril models at pH 7 and pH 2 available, it becomes possible to test the contribution of various charged groups to ΔΔGel. We performed such tests by threading residues with charges appropriate for pH 2, neutral for glutamic acid and +1 for histidines, onto conformations recorded in our simulations at pH 7 and 2. We refer to the residues with altered charge state as “mutated” in this charge scan analysis, although no actual mutation is taking place. For the E22 mutant, conformations obtained at pH 2 were used while for all others - conformations generated at pH 7 with the E22-K16 salt bridge intact. We note that our procedure is far from rigorous and at best can give only an indication as to the effect of a particular mutation. Its predictions need to be tested in actual simulations. Table 2 shows our estimate of the effect that charge alteration of each particular group has on ΔΔGel. Remarkably, the results indicate that only one group out of six, the charge at C-terminus, would reverse the sign of ΔΔGel thus making the registry R2 structure more stable at pH 7. To find a mechanistic explanation for this effect, we probed the interactions of C-terminal charge in the R1 structure at pH 7. Contacts with the side chains of other residues showed a particularly strong interaction between C- and N-terminal groups belonging to two adjacent β-strands. There is only one such interaction per β-strand pair and it occurs alternatively along the β-sheet. Our previous analysis on the hydrogen bond network suggested that this interaction may exist. Here we find that it may play an important role in defining the registry of the fibril. Since β-strands in registry 2 are shifted by 2 residues along each other relative to registry 1, they do not allow for a contact between C- and N-terminal groups to be formed. A cartoon explaining the geometrical difference between R1 and R2 β-sheets is shown in Figure 13. It is seen that in two competing structures R1 and R2, the presence of a C-N contact will shift the balance toward the former. Canceling the C-N interaction makes R2 conformation electrostatically more favorable. This is what Table 2 suggests based on the hypothetical threading experiments. Would this theory be supported by real simulations?
To test the effect of neutralizing C-terminal charge on ΔG, two additional simulations were designed. Conformations R1 – C4– D0 and R2– C4– D4 were prepared in which a neutralizing amino group NH2 was appended to the C-terminal. The simulations were performed for 8ns using the same protocol as described previously. The resulting free energy difference is presented in Table 3, together with the data obtained for the wild type sequence. Comparison of the two data sets shows that the major effect of neutralizing the C-terminus is on the electrostatics free energy component. Although ΔΔGel does not change its sign, as the charge scanning test suggested, its magnitude is reduced by more than 3-fold. With the slightly more favorable hydrophobic interactions in the mutated sequence compared to the original peptide, −106kCal/mol versus −72kCal/mol, the total free energy difference ΔG becomes negative signaling that R2 structure is more stable. Our design hypothesis about the critical role of the C-terminal charge thus seems to be correct: the presence of this charge does seem to be tilting the structural balance toward registry R1 at pH 7. Based on these results, we predict that Aβ11–25 peptide with neutralized C-terminus will assemble into registry 2 fibrils at pH 7. This prediction is amenable to experimental verification.
Although amyloid fibrils differ significantly in their size and microscopic appearance, they have several common structural elements28. First, X-ray diffraction experiments show that all fibrils are rich in cross β-structure; they are made of β-sheets in which β-strands are perpendicular to the fibril axis while the hydrogen bonds linking β-strands run parallel to it. The distance between two adjacent β-strands is 4.7–4.8Å. Second, two or more β-sheets coalesce to form rod-like structures of 25–30Å across, known as protofilaments. The distance between β-sheets in a protofilament varies according to the size of the side chains of amino acids that lie at the interface between two sheets. Third, fibrils typically consist of several protofilaments wound up around a common axis. Figure 1(c) shows a schematic representation of the hierarchical organization of amyloid fibrils up to the level of protofilaments.
Amyloid β peptide fragment Aβ11–25, with the amino acid sequence EVHHQKLVFFAEDVG, forms fibrils that conform to the hierarchical model of amyloid. Through electron microscopy and X-ray diffraction8; 30, this peptide was seen to grow into straight, unbranched, ribbon-like filaments with a periodic modulation indicating a twist along the long axis. The filaments were determined to consist mostly of β-sheets in anti-parallel arrangement while constraints derived from the solid-state NMR measurements allowed for the registry of hydrogen bonds between β-strands to be determined. It was found that at pH 7, residues 17+k of one β-strand are aligned with the residues 20-k of the next β-strand thus forming registry 17 k ↔ 20 − k, or R1 as it will be referred to in this work. At pH 2, the β-sheets are still anti-parallel but their registry changes to 17 + k ↔ 22 − k, or R2. Microscopic details of the two registries are explained in Figure 1(a)–(b). R1 and R2 display different patterns of inter-strand hydrogen bonds, in particular the maximum number of permitted bonds. Side chain contacts differ both at the level of a β-sheet as well as in the context of a fibril. There are likely other characteristics in which the two registries differ, such as the exposure of amino acid side chains to solvent. The exact microscopic structure of Aβ11–25 fibrils remains unknown.
We use all-atom representations to model peptides and water in all our simulations. Fixed-charge OPLS/AA force field31 was applied to model interactions among peptide atoms and TIP3P32 potential was used to model water. All simulations were performed using the GROMACS software set33; 34. Covalent bonds of the water molecules were held constant by SETTLE algorithm35. The bonds involving hydrogen of the peptide were constrained according to LINCS protocol36. The integration time step was 2 fs. Nonbonded Lennard–Jones interactions were tapered starting at 7 Å and extending to 8 Å cut-off. Neighbor lists for the nonbonded interactions were updated every 10 simulation steps.
To compute relative free energy of forming a β-strand dimer in registries R1 and R2, we use the umbrella sampling method18, combined with the multiple-histogram analysis tool19; 20. Initial conformations of the peptides were constructed in the ideal anti-parallel β-sheet configuration with /ψ angles at −119 and 113 values, respectively37. Eight Cl ions were added to the system at pH 2 and four Na ions at pH 7 to render the total charge zero. The particle-mesh Ewald method38 was used to compute electrostatic interactions. Nose-Hoover thermostat39 was applied to maintain the temperature in our simulations close to 300K.
The peptides were initially placed in a β-sheet configuration with their Cα lying in the X–Y plane and inter-sheet hydrogen bonds running parallel to Z axis. The coordinates of Cα of one of the two peptides, the lower β-strand, were restrained to lie around their initial values using harmonic potentials with the force constant of 104 kJ/mol/nm2. Harmonic restraints along X- and Y-directions were also applied to Cα of the second peptide, top β-strand, restricting the allowed motion to Z-axis only. Additionally, an angle restraint was applied to the top β-strand, designed to keep the vector formed by Cα’s of E11 and G25 approximately perpendicular to the vertical axis. The resulting backbone conformations of the top β-strand were parallel to the lower β-strand. Identical sets of constraints were applied at neutral and acidic pH to make sure that the results are consistent under both sets of conditions. Umbrella sampling potentials were defined in terms of the distance between centers of mass of Cα atoms in the two β-strands. Forces resulting from these potentials were equally distributed among all Cα atoms. 19 different windows along Z-axis, starting at 4Å and spaced by 0.5Å, were considered with the force constant of 2×104 kJ/mol/nm2. Umbrella sampling simulations were run for 5ns. The resulting histograms were tested for sufficient overlap before computing potentials of mean force.
We used the same protocol as in our free energy calculations to perform short simulations of the proposed fibril models. The initial conformations were prepared as follows. First, peptides were set up in the ideal anti-parallel β-sheet configuration. Next, they were translated by 4.8Å to form a β-sheet of eight (or 16) peptides. The second β-sheet was copied from the first β-sheet and then translated by 8Å and rotated as necessary. The models were solvated in rectangular boxes of water making sure that no water molecules are trapped in their interior. Simulation boxes were chosen to allow for at least 10Å of buffer layer of water between peptides and the sides of the box. The dimensions of the boxes varied for each model, with a 9×6×6nm box being representative for 16-chain models. Models with 32 peptides were approximately twice as long. The initial configurations were equilibrated for 20ps with the Cα atoms fixed. The equilibration phase was followed by multiple unconstrained production runs of 4ns. In addition to temperature, pressure was also kept constant in these simulations at P=1atm. Electrostatic interactions were treated using the reaction-field method40; 41. Atom-based cutoffs were applied to water molecules and charge group-based cutoffs to peptides to better model hydration of charged groups in this method42. The cut-off distance of 12 Å was used in all simulations.
To differentiate among many fibril models generated in this study we estimate their relative free energy. We start with the general expression for free energy of a protein in state Γ:
where U(Γ) is the internal energy of the protein in vacuum, T is temperature, S(Γ) is the entropy associated with the state and ΔG(Γ) is the cost of solvating the protein in water, or solvation free energy. State Γ in these notations refers to distributions of protein conformations, rather than single conformations. Two components in eq (1) can be further decomposed into more elementary terms: U(Γ)= Uvdw (Γ) + Uel (Γ), where Uvdw (Γ) and Uel (Γ) are the van der Waals and electrostatic interactions among protein atoms respectively, ΓG(Γ) = Gpol (Γ) + ΔGnp (Γ), where ΔGpol (Γ) is the electrostatic, or polar, solvation free energy and ΔGnp (Γ) is the non-polar, or hydrophobic, solvation energy. In this work, we ignore the entropic term in expression (1). Although the magnitude of S(Γ) is certainly non-negligible, we expect it to be similar for similar fibril models that differ only in their registry, and therefore to drop from (1) when free energy difference is taken. We will also group two electrostatic terms ΔGel (Γ) = ΔGpol (Γ) + Uel (Γ) to arrive at the final expression for free energy:
The polar solvation free energy is evaluated in this work using the molecular volume formulation of the generalized Born model (GBMV)43. The GB model is one of the most accurate solvation schemes currently available. It was seen to generate results in very good agreement with the solution of the more accurate Poisson equation44 for a large variety of protein structures. Our present understanding of the non-polar solvation is much less developed than that of the polar solvation. Hydrophobic interactions exhibit a subtle dependence on the size of the solute45 and therefore are not trivial to model theoretically. For sufficiently large solutes, starting at approximately 10Å across, the cost of burying a patch of hydrophobic surface on a solute molecule is proportional to its surface area (SA). β-sheets involved in our fibril models reach the dimensions where the SA model applies. We therefore estimate the hydrophobic cost associated with forming a fibril model as proportional to the surface area ΔS that becomes shielded from solvent when β-sheets are stacked together, ΔGnp (Γ) = γΔS. The proportionality constant is set to the value that corresponds to flat surfaces, γ =72 Cal/mol/Å2 46. All calculations are performed using CHARMM software suite47.
We gratefully acknowledge the support of the North Carolina Biotechnology Center and the National Institutes of Health (to AB), grant # R01GM083600-02
Publisher's Disclaimer: This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.