Molecular dynamics simulation on the effect of transition metal binding to the N-terminal fragment of amyloid-β

Abstract

We report molecular dynamics simulations of three possible adducts of Fe(II) to the N-terminal 1–16 fragments of the amyloid-β peptide, along with analogous simulations of Cu(II) and Zn(II) adducts. We find that multiple simulations from different starting points reach pseudo-equilibration within 100–300 ns, leading to over 900 ns of equilibrated trajectory data for each system. The specifics of the coordination modes for Fe(II) have only a weak effect on peptide secondary and tertiary structures, and we therefore compare one of these with analogous models of Cu(II) and Zn(II) complexes. All share broadly similar structural features, with mixture of coil, turn and bend in the N-terminal region and helical structure for residues 11–16. Within this overall pattern, subtle effects due to changes in metal are evident: Fe(II) complexes are more compact and are more likely to occupy bridge and ribbon regions of Ramachandran maps, while Cu(II) coordination leads to greater occupancy of the poly-proline region. Analysis of representative clusters in terms of molecular mechanics energy and atoms-in-molecules properties indicates similarity of four-coordinate Cu and Zn complexes, compared to five-coordinate Fe complex that exhibits lower stability and weaker metal–ligand bonding.

Communicated by Ramaswamy H. Sarma

Keywords:

• Amyloid-&beta
• peptide
• transition metals
• molecular dynamics
• AMBER
• Ramachandran map

Introduction

Neurodegenerative diseases present one of the greatest health-care challenges to modern society. Alzheimer’s disease (AD) affects millions and is projected to increase markedly in the coming decades (Kepp, 2017; WHO, 2017; World Alzheimer Report 2015, 2015). AD is characterised by degeneration and death of nerve cells and their connections in addition to the presence of insoluble plaques and fibrils in diseased brains (Hardy & Higgins, 1992; Karantzoulis & Galvin, 2011; Kepp, 2012; Kosik, Joachim, & Selkoe, 1986). The amyloid (Cleary et al., 2005; Haass & De Strooper, 1999; Selkoe, 1991; Walsh et al., 2002) and metal ion (Bush & Tanzi, 2008; Kepp, 2017) hypotheses are leading theories to explain the onset of AD: in the former, accumulation of the β-amyloid (Aβ) peptide into toxic soluble oligomers as well as observed plaques is argued to be the cause of disease. In the latter, dysregulation of naturally occurring metals, notably copper, zinc and iron through genetic, environmental or age-related factors, is believed to promote aggregation and hence disease onset (Bush, 2013; Faller & Hureau, 2012; Kepp, 2017; Kozlowski, Luczkowski, Remelli, & Valensin, 2012; Lovell, Robertson, Teesdale, Campbell, & Markesbery, 1998; Pithadia & Lim, 2012); the exact role and behaviour of these metals in relation to AD is reviewed in detail elsewhere (Duce & Bush, 2010; Greenough, Camakaris, & Bush, 2013; Kepp, 2012, 2017; Kozlowski et al., 2012).

Metal ion coordination has distinct effects on the structure and chemistry of Aβ, including its propensity to aggregate. The anchoring effect of metal coordination provides structural restraint and potentially a seed for aggregation. Indeed, many studies report that metal ions act as strong inducers of Aβ aggregation (Bush et al., 1994; Esler et al., 1996; Huang et al., 2004) though the results are diverse (Miller, Ma, & Nussinov, 2012), including a range of fibrillar and non-fibrillar, amorphous aggregates. Furthermore, these results are strongly dependent on reaction conditions. The redox activity of Cu and Fe may also lead to the formation of reactive oxygen species that place diseased nerve cells and synapses under oxidative stress (Greenough et al., 2013).

Aβ contains a high-affinity metal-binding region at residues 1–16, which includes three histidine residues (His6, His13 and His14) as well as other groups (e.g. Asp and Glu) that can coordinate metal ions depending on conditions (Curtain et al., 2001; Faller & Hureau, 2009; Faller, Hureau, & La Penna, 2014; Kepp, 2012, 2017; Kozin, Zirah, Rebuffat, Hui Bon Hoa, & Debey, 2001). These metal–Aβ complexes usually form at 1:1 ratio, though secondary, lower-affinity binding has been reported (Atwood et al., 2000; Kepp, 2017; Kozin et al., 2001; Tõugu, Karafin, & Palumaa, 2008; Tõugu & Palumaa, 2012). To date, most studies have focussed on Cu(II) and Zn(II) binding to Aβ. In contrast, few studies have identified Fe(II)–Aβ interactions – suggested to occur via Asp1, Glu3, His6, His13 and His14 (Bousejra-ElGarah, Bijani, Coppel, Faller, & Hureau, 2011). However, given the seemingly major role of iron dysfunction in AD (Bush, 2013; Honda, Casadesus, Petersen, Perry, & Smith, 2004; Kepp, 2012, 2017), further investigation of iron–Aβ complexes is required.

Molecular simulation is a well-established tool to examine metal interactions with Aβ: many studies using a range of methods, from ab initio and DFT through QM/MM to atomistic molecular dynamics (MD), have been reported for both Cu(II) and Zn(II) (Alí-Torres, Maréchal, Rodríguez-Santiago, & Sodupe, 2011; Alí-Torres, Mirats, Maréchal, Rodríguez-Santiago, & Sodupe, 2015; Furlan, Hureau, Faller, & La Penna, 2010; Giannozzi et al., 2012; Huy, Vuong, La Penna, Faller, & Li, 2016; La Penna, Hureau, Andreussi, & Faller, 2013; Pan & Patterson, 2013; Raffa & Rauk, 2007; Wise-Scira, Xu, Perry, & Coskuner, 2012). However, comparisons of the effects of different metals on Aβ under comparable simulation conditions are uncommon; to the best of our knowledge, only one report of MD simulations on Fe(II)–Aβ systems has been reported in the literature (Boopathi & Kolandaivel, 2016), so our initial focus is on Fe(II) bound in three different proposed coordination modes. These are shown in Figure 1: Fe_1, where Fe(II) is bound via Asp1 (NH₂ and CO), His6 (Nδ and CO) and His14 (Nδ); Fe_2, where Fe(II) is bound via Asp1 (NH₂ and CO), His6 (Nδ and CO) and His13 (Nδ); and Fe_3, where Fe(II) is bound via Asp1 (NH₂ and CO), Glu3 (Oε), His6 (Nδ and CO) and His14 (Nδ). Then, to compare the binding of different biologically relevant transition metals to Aβ, we have carried out microsecond timescale MD of the N-terminal 1–16 fragments of Aβ to Cu(II), Fe(II) and Zn(II), taking binding modes from those proposed by Nasica-Labouze et al. (2015) – illustrated in Figure 1 – specifically, Cu(II) bound through Asp1 (NH₂ and CO), His6(Nε) and His14(Nε); Fe(II) bound through Asp1 (NH₂ and CO), His6 (Nε and CO) and His14 (Nε); and Zn(II) bound through Asp1 (Oδ), His6 (Nε), Glu11 (Oε) and His14 (Nε). The purpose of the current study is twofold: firstly, to examine whether changes in the identity of the particular amino acid residues bound to the metal strongly affect the peptide structure obtained from MD, using three different binding modes of Fe(II) as a test case. We then compare the effect of different metals on peptide structure and dynamics using comparable binding modes. Our primary focus in this work is on establishing a suitable protocol for MD simulation of these flexible systems: we hope to report simulations of other important systems, such as trivalent Fe(III) or Al(III) ions or non-natural metals such as Pt(II) or Co(III), in future publications.

Figure 1. Metal binding modes used in this study. Top row, left-right: Fe_1, Fe_2 and Fe_3. Bottom row, left-right: Cu and Zn.

Computational methods

Aβ1–16 was constructed in an extended conformation in MOE (Molecular Operating Environment (MOE) (2013.08), 2013) with the appropriate residue protonation states for physiological pH. Metals were coordinated to the peptide as shown in Figure 1 and subjected to brief LowMode conformational searches (Labute, 2010) to obtain reasonable starting structures. MD simulations were performed using the AMBER16 package (Case et al., 2016). The AMBER ff14SB forcefield parameter set (Maier et al., 2015) was used to model all standard amino acid residues, while metal parameters were obtained using the MCPB.py program (Li & Merz, 2016). Here, parameters are obtained from B3LYP/6-31G(d), and RESP charges for the metal-coordinating regions were obtained at the same level of theory using Gaussian09 (Frisch et al., 2009). The geometry of each system was optimised using 1000 steps of steepest descent and 1000 steps of conjugate gradient methods. MD simulations were carried out in the NVT ensemble, using a Langevin thermostat to control the temperature at 310 K. The generalised Born solvation model (Clark Still, Tempczyk, Hawley, & Hendrickson, 1990; Constancien & Contreras, 1984; Schaefer & Karplus, 1996) was used to solvate all systems, with electrostatic interactions neglected beyond a cut-off of 12 Å: this approach has been shown to enhance conformational sampling of flexible systems (Anandakrishnan, Drozdetski, Walker, & Onufriev, 2015). During all simulations, the SHAKE algorithm (Ryckaert, Ciccotti, & Berendsen, 1977) was used to constrain bonds to hydrogen. Simulations were performed using a 2-fs integration time step. Analysis of the trajectories was performed using CPPTRAJ v16.16 (Roe & Cheatham, 2013) and VMD 1.9.3 (Humphrey, Dalke, & Schulten, 1996). Secondary structure assignment was made using the DSSP algorithm (Kabsch & Sander, 1983) via CPPTRAJ. Salt bridges were defined as any contact distance of less than 3.2 Å between oxygen and nitrogen atoms in charged residues. Ramachandran maps and hydrogen-bonding plots were made using MDplot (Margreitter & Oostenbrink, 2017). Atoms-in-molecules (AIM) analysis used the AIMAll package (Keith, 2017).

3. Results and discussion

3.1. Effect of Fe(II) coordination mode

Three separate 500 ns MD simulations of each Fe–Aβ complex (denoted Fe_1, Fe_2 and Fe_3) were carried out, starting from the same minimised structure but with different initial velocities, randomly sampled from the Maxwell–Boltzmann distribution at 310 K. Root-mean-square displacement (RMSD) relative to starting structure was used as our primary measure of equilibration (Huy et al., 2016; Mutter, Turner, Deeth, & Platts, 2018; Turner, Mutter, Deeth, & Platts, 2018). Plots of backbone RMSD for each run are shown in Figure 2, which show that simulations reach stable values after 100–300 ns; equilibration points for each simulation are included in Table S1, supplementary material. All analyses reported are taken from data extracted from frames after each equilibration point, averaged over three separate runs, leading to 900, 950 and 950 ns of simulation data for Fe_1, Fe_2 and Fe_3, respectively. ‘Spikes’ in RMSD are seen in some plots, for instance around 350 ns in Fe_1C: closer inspection indicates that these correspond to relatively minor changes in backbone conformation for a few residues (Figure S1, supplementary material), without affecting metal coordination or the overall structure of the peptide.

Figure 2. Backbone RMSD values (Å) for the simulations of Fe_1, Fe_2 and Fe_3.

Table 1. Total percentage of secondary structure, by type, for Fe(II)–Aβ complexes.

Secondary structure	Helix	Sheet	Other
Fe_1	15.49	0.23	84.28
Fe_2	17.20	0.09	82.71
Fe_3	18.52	0.04	81.44

Averages and standard deviations of RMSD, relative to a common starting point, over these trajectories (Table S2, supplementary material) further confirm equilibration: averages are in the range 2–4 Å with small standard deviations 0.2–0.45 Å. Radius of gyration (R_g) values are similar, with small standard deviations of 0.15–0.2 Å (Table S2, supplementary material). Averages for the three coordination modes converge to similar values of between 7.2 and 7.6 Å, with the largest value observed for the six-coordinate Fe_3. Root-mean-square fluctuations (RMSF) for individual amino acids exhibit similar behaviour, as shown in Figure 3 (data in Table S3, supplementary material). Metal-binding residues exhibit the lowest RMSF values – the most mobile being Glu3 in Fe_3, with an RMSF of 5.69 Å – while larger values of up to 10 Å are found for aromatic residues Phe4 and Tyr10, and also for the C-terminal Lys16 residue.

Figure 3. RMSF data (Å) for three different Fe coordination modes.

Table 2. Percentage of secondary structure types observed over equilibrated trajectories.

SS %	Helix	Sheet	Other
Cu	16.28	0.00	83.72
Fe_1	15.49	0.23	84.28
Zn	22.19	0.00	77.81

Table 3. Statistics of hydrogen bond occupation.

# Hydrogen bonds	Ave	SD	Min	Max
Cu	5.32	2.17	0	15
Fe_1	7.33	1.90	0	15
Zn	4.55	1.65	0	9

Metal coordination can exert profound effects on the structure of Aβ, which has an important role in its neurobiology. Figure 4 shows secondary structure by residue for each Fe(II)–Aβ system, as well as Ramachandran maps for the equilibrated trajectories; Table 1 shows total peptide secondary structure percentage. Plots are very similar, exhibiting large amounts of bend structure in residues 3–11 and a mix of α- and 3₁₀-helix for residues 12–16. The six-coordinate structure favours turn over bend in the N-terminus, but is otherwise similar. Very small amounts (<0.25%) of β-strand structure are found to be present in all three systems. Ramachandran maps reveal that all three complexes have the greatest concentration of conformations in the region expected for α- and 3₁₀-helix, along with smaller but significant populations of β-strand and left-handed helix, the former being more prevalent for the two five-coordinate complexes (Fe_1 and Fe_2). Residue contact maps are also broadly similar, exhibiting relatively short contacts between coordinated residues and longer distances between central and C-terminal ones.

Figure 4. Secondary structure plots, Ramachandran and residue contact maps for Fe(II)–Aβ complexes.

Clustering analysis of all equilibrated trajectories on the basis of peptide backbone dihedrals was also performed, leading to the data in Figure 5, including the minimised structures of the most representative clusters. For both Fe_1 and Fe_3, a single cluster accounts for over 90% of conformations, while for Fe_2, the dominant cluster accounts for 57% of simulation. All three show similar properties, with coil/bend/turn throughout the N-terminal region and a small amount of helix towards the C-terminus. These aspects of the secondary structure are summarised in Table 1, which reports that each simulation spends between 15 and 20% of the time in helix-like conformations, less than 0.25% in strand-like, and over 80% in turn/bend/coil conformations. Overall, we surmise that the specifics of coordinating residues do not strongly affect the conformations visited by or the flexibility of the peptide chain.

Figure 5. Minimised geometry and clustering percentage of three Fe(II) complexes.

3.2. Comparison of Cu, Fe and Zn binding

We therefore decided to compare data from the Fe_1 coordination mode to analogous structures bound to Cu(II) and Zn(II), since these share many of the same binding residues. All three metals are bound through Asp1, His6 and His14, with Glu11 as an additional binding residue for Zn(II) only. Figure 6 shows backbone RMSD values for each simulation, while Table S4, supplementary material, reports RMSD and R_g values for post-equilibration frames for the three systems studied. In each case, RMSD values are relatively small (2–3.5 Å) with small standard deviations (0.2–0.35 Å). Maximum values are very slightly larger for Zn than for Cu and Fe, but differences are rather small. A similar pattern is seen for R_g data, in which Fe is notably more compact than the others and Zn slightly larger than Cu. These data suggest that the five-coordinate nature of the Fe complex restricts the peptide to a smaller, more compact set of structures. As expected, given the harmonic bonding model used for metal coordination, bond lengths involving metal ions oscillate around equilibrium values of around 2.0 Å with standard deviations of ca. 0.05 Å. Bond angles exhibit a similar pattern, exhibiting standard deviations of 2–3°, and reflect the approximately square-pyramidal, distorted square-planar and tetrahedral coordination geometries of Fe, Cu and Zn, respectively. Full data are reported in Table S5, supplementary material.

Figure 6. Backbone RMSD values (Å) for the simulations of Cu, Fe_1 and Zn.

Table 4. Summary of bond critical point propertiesa for metal–ligand bonds (au).

Metal	Ligand	ρ	∇^2ρ	ε	G	V	H
Cu	N (Asp1)	0.085	0.249	0.013	0.097	−0.131	−0.035
	O (Asp1)	0.085	0.363	0.030	0.118	−0.146	−0.028
	N (His6)	0.098	0.339	0.047	0.126	−0.168	−0.042
	N (His14)	0.094	0.316	0.063	0.118	−0.158	−0.040
Fe	N (Asp1)	0.079	0.419	0.371	0.113	−0.122	−0.009
	O (Asp1)	0.075	0.432	0.893	0.112	−0.115	−0.004
	O (His6)	0.066	0.421	1.165	0.105	−0.106	0.000
	N (His6)	0.085	0.501	0.702	0.134	−0.143	−0.009
	N (His14)	0.079	0.482	0.483	0.125	−0.130	−0.005
Zn	O (Asp1)	0.074	0.293	0.047	0.096	−0.119	−0.023
	N (His 6)	0.075	0.236	0.068	0.087	−0.114	−0.028
	O (Glu11)	0.053	0.205	0.063	0.061	−0.070	−0.010
	N (His14)	0.082	0.256	0.057	0.098	−0.131	−0.034

^aρ is the electron density at the bond critical point; ▿2ρ is the Laplacian of the electron density here; ε is the ellipticity of the bond; G is the kinetic energy density at the bond critical point; V is the potential energy density here; H is the total energy density, the sum of G + V.

All properties are reported in atomic units.

Further information on the impact of metal coordination on peptide conformational freedom can be drawn from RMSF data, as reported in Table S6, supplementary material, and Figure 7. The overall pattern for all three metals is similar: metal-bound residues are the least mobile: His6 and His14 in particular exhibit lower RMSF values. Reflecting the values above, Fe_1 has the lowest average RMSF and the lowest value for 9 of the 16 residues, while Zn has the highest average RMSF and is highest for nine residues.

Figure 7. RMSF values (Å) for Cu, Fe_1 and Zn structures.

Metal binding has important effects on the structure and chemistry of Aβ: Figure 8 reports a summary of the secondary structure across equilibrated trajectories, and Table 2 shows total peptide secondary structure percentage. Secondary structures show a good deal of similarity between metals: amino acids Asp1–Tyr10 exhibit mainly coil/turn/bend structures, while Glu11–Gln15 show a much more helical character, along with turn. Within this broad picture, however, subtle details differ: Zn induces more α-helical content than the other metals, which extends further towards the C-terminus of the peptide, whereas Cu and Fe show more 3₁₀ helix characters over a slightly shorter range. A small amount of parallel β-strand contact between Ala2 and Gln15 is also seen for Fe_1, and an even smaller extent of anti-parallel strand for Ser8 and Gly9 in the Cu complex. Ramachandran maps reflect this pattern, but also show more detail of conformations sampled by MD. All show the greatest concentration within the right-handed helix region of the maps, along with variable amounts of β-strand and left-handed helical conformations. Notably, Cu shows a large number of conformations in the PPII region (ca. −75, 150) (Adzhubei & Sternberg, 1993; Hollingsworth & Karplus, 2010; Woody, 2009). The PPII structure is relatively open and contains little/no internal hydrogen bonding (vide infra) and as such is likely characterised as turn/bend/coil. Fe_1 and Zn show progressively fewer conformations in this region. In addition, Fe_1 shows significant occupation of the bridge and ribbon regions of the plot, extending between the broad helical and β-sheet regions (Hollingsworth & Karplus, 2010; Karplus, 1996) Fe_1 also has a greater number of conformations in the β-sheet region of the Ramachandran map than either Cu or Zn, though this is not reflected in secondary structure analysis, likely due to their lack of requisite hydrogen bonds (vide infra). Overall though, we stress that the N-terminus of Aβ remains rather unstructured despite the metal coordination imposing structural restraints on certain residues, as exemplified by the data in Figure 8.

Figure 8. Breakdown of secondary structure and Ramachandran maps for Cu, Fe_1 and Zn complexes with Aβ1–16.

Salt bridges are crucial in determining peptide structure and stability: Figure 9 shows the percentage occupancy of all eight possible such contacts within Aβ_1–16 (numerical data in Tables S7–S9, supplementary material). In all systems, the most commonly observed interaction is between Arg5 and Asp7, despite the fact that the metal ions are coordinated through His6 between these residues. This bridge is present for almost 100% of the sampled configurations of the Fe_1 complex, and between 75–85% for Cu and Zn. Glu3–Lys16 is also very prevalent (91%) for Fe_1, but is barely observed for Cu and Zn for which Glu3 tends to make contact with Arg5 (38 and 10%, respectively). Notably, the Asp7–Lys16 interaction is observed much more frequently (14%) in the Cu simulations than in either the iron or zinc. Asp1 and Glu11 rarely form salt bridge interactions in any of the metal systems.

Figure 9. Occupancy of salt bridges and contact maps for Cu, Fe_1 and Zn complexes.

The effect of these salt bridges, and of the metal coordination patterns, is evident in the contact maps shown in Figure 9. Fe_1 contains rather short distances for distant amino acids that share coordination, for instance, between Asp1 and His14 and also between Glu3 and Lys16 that participate in a high-occupancy salt bridge. In contrast, Cu and Zn complexes have rather fewer close contacts, reflecting the lower coordination number and fewer salt bridges in these complexes. Particularly long distances are observed for Glu3–Val12 in Cu, and Glu3–Tyr10, His13 and Lys16 in Zn.

As well as salt bridges, hydrogen bonds play an important role in determining stability and structure: Table 3 summarises the overall H-bond patterns for the three metals. All are rather transient, with frames containing zero hydrogen bond interactions observed as well as others with as many as 15. On average, Fe_1 contains the most: contacts between backbone C = O of Glu11 and H-Nε of His14 (58% of trajectory) and between backbone C = O of Phe4 and H-Nε of His6 (57%) are the most prevalent. In contrast, the most occupied H-bond in the Zn complex is between Tyr10 C = O and H-Nε of His14 at just 33% of frames, while for Cu, the contact between Lys16 C = O and H-Nε of His6 is observed for 42% of frames. Figure 10 adds more detail to this picture and highlights significant differences between metals: both Cu and Zn contain frequent interactions between Arg5 donor and Glu3 and Asp7, whereas Fe contains only the former interaction. Both Cu and Fe display frequent His6 (donor) to Glu3 and Phe4 (acceptors), which are not present in Zn. In contrast, Zn Glu11 accepts H-bonds with four different residues with greater than 5% incidence, the most of any residue. Interestingly, Glu11 is coordinated to Zn in this case, while the uncoordinated Glu11 in Fe and Cu cases forms H-bond interactions with no more than three residues, though at greater incidence than for Zn.

Figure 10. Plots of hydrogen bond occupancy for Cu (left), Fe_1 (centre) and Zn (right). Plots only display data for hydrogen bonds present for greater than 5% of simulation. Black circles indicate multiple interactions between relevant residues.

Clustering of all recorded frames on the basis of backbone dihedrals was also performed, leading to the data shown in Figure 11, including ribbon diagrams of the most representative cluster following minimisation. For Cu and Zn, a single cluster accounts for over 60% of the trajectory, with the remaining frames spread across clusters encompassing no more than 17%. The Fe trajectory, however, leads to essentially a single cluster accounting for the entire set of post-equilibration frames recorded at 92%, with no other cluster covering more than 5%. This reinforces the picture that the extra metal–ligand bond in the five-coordinate iron complex limits the flexibility of the peptide more than is the case for four-coordinate metals. We also report the total MM energy after minimisation, which shows that the Cu and Zn complexes have similar energy and are markedly more stable than the Fe_1 complex. Further inspection of the breakdown of this energy (Table S10, supplementary material) indicates that bond stretch, angle bend and van der Waals terms are within 2 kcal/mol for all three systems and that the overall energy difference stems from changes in electrostatic and generalised Born solvation terms. It seems that higher coordination number of the iron complex prevents the polar side chains from interacting electrostatically with either the rest of the peptide or the implicit solvent, compared to the four-coordinate Cu and Zn complexes.

Figure 11. Summary of clustering analysis: minimised geometry, percentage occupancy and molecular mechanics energy (kcal/mol) after minimisation.

Optimised geometries of representative clusters were used to carry out AIM analysis, based on electron density from B3LYP/6-31G(d) calculations. This is focussed on the metal–ligand bonds: bond critical points (BCPs) are located for all coordination bonds: properties at these points are summarised in Table 4. This reveals a pattern of metal–His bonds being notably stronger than those of other atoms, exhibiting greater electron density (ρ), lesser positive Laplacian (∇²ρ) and more negative energy density (H). By these measures, the strongest of these bonds are between Cu and Nδ of His6 and His14, while those to Asp1 are around 10% weaker. This pattern is less marked for the other metals, but in general, these data indicate that coordination through oxygen is weaker than through N. Particularly weak, and hence more labile, bonds are observed for Fe—O(His6) and Zn—O(Glu11).

4. Conclusions

We have used atomistic MD to examine and compare the effect of coordination of copper, iron and zinc on the structure and dynamics of the N-terminal 1–16 fragments of Aβ. These ions are believed to play a role in the aggregation of Aβ and in the onset of AD, such that their effect on peptide conformational behaviour may be an important factor in disease. In order to compare different metals and coordination modes on an equal footing, we use Li and Merz’s MCPB.py approach, which extracts force constants and non-bonded parameters for metal ions from DFT calculations. These parameters were calculated for three different proposed coordination modes for Fe(II) and for comparable coordination of Cu(II) and Zn(II).

This method is used to obtain 1.5 μs of MD trajectories for two five-coordinate Fe(II) complexes, which differ in coordination through either His13 or His14, and a six-coordinate complex involving also the side chain carboxylate of Glu3. Each system is found to reach pseudo-equilibration after 100–300 ns, allowing us to extract over 900 ns of equilibrated data. Over the course of these trajectories, the conformations adopted by the peptide are broadly independent of the coordination mode: all exhibit turn/coil/bend structures in the N-terminus and a significant amount of helical content towards the C-terminus. Similar patterns in the flexibility of individual residues, as measured by RMSF, and inter-residue contact distances are also seen.

Comparison of one of these Fe(II) coordination modes with analogous ones to Cu(II) and Zn(II) was then made, from which we find more notable differences between metals. Equilibration again occurs within similar time scales, allowing extraction of at least 900 ns of equilibrated trajectories. The overall picture of turn/bend/coil in N-terminus and helix in C-terminus is retained, but the amounts of each vary more between metals. In particular, Cu coordination leads to greater occupancy of the PPII region of Ramachandran maps more than other metals, while Fe leads to significant occupation of bridge and ribbon regions. Salt bridges are also affected differently by different metals: Arg5–Asp7 is ubiquitous, whereas Glu3–Lys16 is common for Fe complexes, but is hardly present for Cu/Zn. These changes are evident in inter-residue contacts: the five-coordinate Fe complex is more compact than the four-coordinate Cu/Zn ones, with Zn in particular allowing greater peptide flexibility. Potential biological consequences of these changes are difficult to estimate, especially since we have only considered the N-terminal portion of the peptide. Nevertheless, we speculate that the rigidity of the Fe complex might provide a better seed for aggregation than the more flexible Cu and Zn complexes.

Clustering allows us to extract a single structure that is representative of the equilibrated trajectories: after minimisation, the energies of the Cu and Zn complexes are similar, but that of the Fe complex is significantly higher, apparently due to the strain imposed by the extra coordination bond for this metal and the lower solvation energy that results from coordination of the carboxylate of Glu3. DFT and AIM analysis of the metal–peptide bonds within these minimised clusters indicates that metal–histidine bonds are strongest, especially within the Cu complex, and that bonds to N are slightly stronger than those to O.

This work focusses on the N-terminal fragment of Aβ since this is known to be the site of metal binding, while the relatively compact nature of the complexes formed aid equilibration of MD trajectories. Having established suitable parameters for simulation of different transition metals, future work will address how the same methods can be applied to the full 1–40 or 1–42 sequences and their aggregation.

Data availability statement

The data that support the findings of this study are openly available via Zenodo at 10.5281/zenodo.1540579.

Acknowledgements

We are grateful to Cardiff University’s ARCCA facility for access to high-performance computing.

Disclosure statement

No potential conflict of interest was reported by the authors.

Additional information

Funding

This work was funded by the UK Engineering and Physical Sciences Research Council (EPSRC) under grant reference [EP/N016858/1].