Structural and Spectroscopic Properties of Assemblies of Self-Replicating Peptide Macrocycles.

Self-replication at the molecular level is often seen as essential to the early origins of life. Recently a mechanism of self-replication has been discovered in which replicator self-assembly drives the process. We have studied one of the examples of such self-assembling self-replicating molecules to a high level of structural detail using a combination of computational and spectroscopic techniques. Molecular Dynamics simulations of self-assembled stacks of peptide-derived replicators provide insights into the structural characteristics of the system and serve as the basis for semiempirical calculations of the UV-vis, circular dichroism (CD) and infrared (IR) absorption spectra that reflect the chiral organization and peptide secondary structure of the stacks. Two proposed structural models are tested by comparing calculated spectra to experimental data from electron microscopy, CD and IR spectroscopy, resulting in a better insight into the specific supramolecular interactions that lead to self-replication. Specifically, we find a cooperative self-assembly process in which β-sheet formation leads to well-organized structures, while also the aromatic core of the macrocycles plays an important role in the stability of the resulting fibers.

Nowadays, living systems replicate using complex protein machinery that copies information-rich molecules such as DNA.[1] However, near the origin of life, simpler molecular self-replication in which molecules make copies of themselves is believed to be a more relevant mechanism.[2−4] Only when exponential growth is achieved, such molecules can outcompete others for common building blocks leading to survival of the fittest and extinction of the weakest replicators.[5−9] One of the main mechanisms through which exponential self-replication has been realized is self-assembly of the replicating species into one-dimensional nanostructures under continuous mechanical agitation.[9−11] While these structures grow from their ends, mechanical energy can break the fibers into two parts, doubling the number of growing ends. The capability of the longer fibers to absorb and dissipate mechanical energy drives the system to a self-assembled state that can be away from thermodynamic equilibrium.[12]

In order to fully comprehend the mechanism of self-replication and guide the design of self-replicating molecules in the future, we need to understand how and why replicating molecules reach a stable self-assembled state, while their precursors do not. A well-studied example of such a system is the functionalized peptide macrocycles described by Otto and co-workers.[9,11,13] We focus on the a pentapeptide (1) with an aromatic dithiol headgroup (Figure ) which is oxidized (either by sodium perborate oxidation or by slow atmospheric oxygen oxidation) to form an interconverting mixture of macrocycles (a “library”) around neutral pH through disulfide exchange reactions. For the peptide sequence 1, it has been found that under constant agitation (stirring of the reaction mixture), the hexamer macrocycle self-replicates driven by assembly into nanofibers (Figure ) and, consequently, the concentration of the hexamer grows exponentially at the expense of other macrocycles (i.e., trimer, tetramer).

Figure 1

Chemical structure of the peptide monomer and schematic representation of the system. Eventually, the collection of different macrocycles (trimer, tetramer, hexamer and other ring sizes) is completely converted to hexamers with self-assembly as the driving force. For more details see, e.g., refs (9, 11).

Several theoretical and experimental methods are able to interrogate the structural and optical properties of the assembly of replicating species. Here, we focus on molecular dynamics (MD) and electron microscopy to study the architecture of the nanostructures, while Fourier transform infrared (FTIR) absorption, UV–vis absorption and circular dichroism spectroscopy are used in combination with semiempirical computational methods to understand peptide secondary structure, aromatic stacking and chirality, respectively.

MD simulations at various levels of detail are routinely used to study the structure and dynamics of biomolecules, such as proteins and lipids. With the recent advances in computational power, simulations of the properties of large self-assembled bioinspired nanostructures are becoming accessible.[14−24] On the other hand, modeling the optical properties of self-assembled nanostructures based on MD trajectories is a developing field that requires different techniques in the different ranges of the light spectrum.[25,26] Considering that the properties that these structures acquire originate from the collective action of many molecules, quantum chemical techniques that model the entire electronic structure of the assembly are very computationally demanding. However, semiempirical methods have been developed to allow such calculations. Specifically for the determination of the secondary structure of proteins and peptides, well-tested frequency mapping algorithms[27−29] have been developed to construct the vibrational Hamiltonian of the system and thus predict the infrared absorption of all the peptide bonds in a molecule. Interestingly, these methods are not typically applied to the self-assembly of short peptides, which make up a large class of biocompatible materials.[30] Because of the many identical molecules in a self-assembled nanostructure, these systems typically represent a well-averaged ensemble and should thus form an ideal case for the frequency mapping algorithms originally developed for proteins. In addition, the structural heterogeneity of supramolecular systems (and sometimes polymorphism[16,31,32]) is intrinsically sampled when multiple snapshots from the MD are considered in higher level calculations.

In the vibrational signature of the peptide-based nanostructures, the amide I mode is the most important secondary structure marker, but in simulations of absorption in the UV and visible range of the spectrum more than one vibronic excitation needs to be considered due to the high density of states in systems with many chromophores. Specifically the absorption of circularly polarized light is of relevance, as individual chromophores might not exhibit a CD signal, but a chiral supramolecular organization often induces such a component. Spano et al. have proposed an INDO (intermediate neglect of differential overlap[33]) based method to calculate UV–vis and such supramolecular CD signals from first-principles (both absorption and emission) and shown its application to 1D aggregates of oligophenylenevinylene.[34−36] As it was previously hypothesized that the peptide macrocycles also form helical assemblies,[13] the comparison of the calculated supramolecular signal with the experimental signal is a good test for the proposed models.

In this paper, we will first discuss the results on the nanoscale architecture of the 1 monomers and fibers made of 1 as observed by MD and transmission electron microscopy (TEM). The resulting structure is then investigated in even greater detail by comparing experimental and computational spectroscopic data that interrogate specific intermolecular interactions.

Results and Discussion

Structural Studies Reveal One-Dimensional Assembly

Libraries of 1 were prepared in two different conditions. In the first case, the initial dithiol monomer is oxidized to 80 mol % using sodium perborate solution. In the second case, the dithiol is left to slowly oxidize by exposure to air. Both libraries were confirmed to contain >95% 1 after 2 weeks. TEM images of the library of 1 indicates the self-assembly of hexamer macrocycles into 1D fibers with a diameter of 3.4 ± 0.7 nm (preoxidized library) and 3.1 ± 0.5 nm (air-oxidized library) and a length ranging up to hundreds of nanometers (see Figure and Figure S1, S2). These measurements suggest the fibers consist of a stack of single macrocycles (cf. Figure C). Interestingly, almost all observed structures in the TEM image of the air-oxidized library are bundles of two individual fibers with a left-handed helical pitch of 56 ± 8 nm, which is similar to matured cross-β amyloid fibrils, which are reported to aggregate by a minimum of two and had a helical pitch of 85–100 nm.[37,38]

Figure 2

Negative stain TEM images of hexamer fibers. (a) Single fibers in a chemically preoxidized sample. Scale bar: 100 nm. Fibers exhibit a diameter of 3.4 ± 0.7 nm. (b) Fibers in bundles of 2 in the slowly oxidized sample, with a width of 3.1 ± 0.5 nm per fiber and a helical pitch of 56 ± 8 nm. Scale bar: 50 nm.

Figure 3

Chemical structure of hexamer macrocycle in (a) the “cartwheel” conformation and (b) the “pairwise” conformation. (c,e) Top and side view of the results of a 100 ns MD simulation of a stack of 16 macrocycles in the cartwheel conformation. The aromatic headgroup is displayed in licorice representation, while the peptide parts are displayed in the ribbon representation using the VMD program. In the side views, all side chains are omitted for clarity, except from the red phenylalanine side chains in (e) to indicate the chiral packing of the stack. (d,f) Top view and side view of the results of a 100 ns MD simulation of a stack of 16 macrocycles in the pairwise conformation. Phenylalanine side chains indicate the hydrophobic “zipper”. Water molecules are not shown.

Typically, aromatic peptide amphiphiles assemble into so-called π–β structures, where the aromatic parts of the molecule stack via π–π interactions and the peptide moieties assemble into parallel or antiparallel β-sheets via main-chain H-bonding and side group interactions.[30] This is likely to be the case for peptide macrocycles, so stacking and H-bonding have to occur along the fiber axis, as will be shown below. Following these constraints, we propose two distinct conformations based on the TEM observations: the C6-symmetric “cartwheel” conformation (Figure A) where all peptide strands are oriented in the same direction, or the C3-symmetric “pairwise” conformation where alternating peptide strands are rotated by 180° to pair the hydrophobic faces (with amino acids Leu and Phe) of the peptide strands (Figure B). Both parallel[39−41] and antiparallel[42−44] conformations are commonly observed among self-assembling short peptides, but since the peptides are connected to each other via an aromatic moiety and the fiber width is equal to the molecular dimensions of the macrocycle, the structures must resemble peptide amphiphiles exhibiting parallel β-sheets such as those designed by Stupp and co-workers.[15,45]

We tested the stability of these conformations by MD simulations and found that both conformations undergo structural changes with respect to their starting structure (see Supporting Information, Figure S3), but the one-dimensional nature of the assemblies is maintained at room temperature (Figure C–F), although partial desorption at the fibril ends can occur. Notably, water is immediately dispelled from the core of the macrocycles leading to the formation of a collapsed fibril with lower symmetry. Lateral sheet–sheet interactions stabilize the structure in both conformations, as can be seen in the Figure C,D. Moreover, the relative packing distances of the aromatic core (parallel displaced stacking with a distance of 3.8–4.2 Å) and the peptide strands (parallel stacking with a Cα–Cα distance of 4.8–4.9 Å) cause the structure to form a helical pitch as visible in Figure E. Interestingly, simulations of 1 consistently formed left-handed helical fibers, while simulations of its enantiomer with d-amino acids consistently formed right-handed helices as can be seen in Figure S4. Figure S5 displays the result of a simulation of two left-handed helical fibers forming a possible conformation of a left-handed bundle, indicating that the correct handedness and pitch (calculated: 46 nm, experimental: 56 nm) for the fibers could be observed in the simulations. However, as there are multiple ways to construct the bundle of fibers and this significantly increases the computational cost of all calculations, we focus on the single fibers in this manuscript for our calculations. All data experimental pertaining to the double fibers of 1 can be found in the Supporting Information.

The cartwheel conformation is stabilized by salt-bridges between the NH3+ groups of Lys3 of the pentapeptides and the C-termini of the neighboring β-sheet (Figure C), while the side chains of Lys5 point into the solution. Leucine and phenylalanine side chains are efficiently shielded from the solvent, except from partial interfaces at the top and bottom of Figure C, which could lead to fiber bundling as observed in the TEM image. For the pairwise conformation, all lysine side chains are solvated and the Leu and Phe side chains are mostly buried in the hydrophobic “zipper” (Figure D,F).

The relative contributions to the stability of the assembly from the aromatic core and the peptides were determined using umbrella sampling simulations as complete binding events were not observed during unbiased simulations. An unbinding trajectory was first generated by pulling a single hexameric macrocycle consisting of only the peptide headgroups truncated at the first amide bond (see Figure S6) from a short stack of either 8 or 2 hexameric headgroups. This trajectory was divided into so-called windows where the centers of mass of the cycle were 0.12 nm apart along the unbinding direction. Sampling the conformational space for 10 ns in the 36 windows thus generated and combining the trajectories using the weighted histogram analysis method (WHAM)[46] resulted in a binding free energy ΔGbind = −70 ± 4 kJ/mol for the stack of 8 macrocycles, and ΔGbind = −50 ± 5 kJ/mol for the stack of two macrocycles, indicating a strong cooperative binding effect. Note that the absolute values of the binding free energies reported here seem to be overestimated, perhaps due to undersampling of the conformational freedom of the macrocycles in solution.

Interestingly, when a similar procedure was attempted for a single β-sheet of 1 with an acetylated N-terminus instead of the dithiol headgroup, it was observed that this sheet is not stable at room temperature and spontaneously dissociates within the 10 ns sampling windows. The binding free energy of the peptides was therefore calculated using simple counting of dimers in a 1 μs simulation of two peptides and found to be −12 ± 0.4 kJ/mol per pentapeptide chain (for details see Supporting Information). The fact that the contribution to the binding free energy of a single peptide is relatively small is not very surprising as the isolated peptide is soluble in water without forming aggregates. Yet, the cooperative assembly of six peptide chains contained in a hexamer can still contribute up to −72 kJ/mol comparable to the stacking interaction of the aromatic core. Moreover, the orientation of the peptide chains appears to be crucial in determining the orientation of the macrocycles as they stack into a one-dimensional fiber.

The MD simulations allow us to speculate on the mechanism of self-replication of one specific size of macrocycle in the experimental library: the cooperative nature of the assembly and the strong negative binding free energy for both peptide and headgroup components imply that nucleation of (especially small) macrocycles into a protofibril is difficult, but the nucleus, once formed, is very stable. This is in good agreement with previous observations that stacking is only observed beyond a certain minimum ring size,[13] which depends on the hydrophobicity of the peptide chain (and therefore on the strength of the peptide–peptide interactions). If we reverse the unbinding trajectory observed in the umbrella sampling simulations, we can conclude that binding of a new macrocycle to the fiber likely occurs via the peptide component first, as the aromatic center of a single macrocycle is shielded from water as much as possible by the peptide arms. Unfolding of the single cycle then leads to rapid attachment. Macrocycles of a different size will not be able to fulfill all favorable hydrophobic or dipole–dipole interactions and can therefore unbind again, thus favoring self-replication under dynamic covalent conditions. Growth via defects in the fiber, as observed in some other systems,[47] is less likely considering the strongly negative binding free energy and a previous demonstration of growth from the fiber ends.[48] We denote that in our unbiased simulations single macrocycles could bind loosely to the side of the fibers, which can act as a mechanism to speed up fiber growth via 1D diffusion along the fiber to the ends, which will be the focus of future work.

Optical Properties Confirm π–β Interactions

As both conformations of fibrils of 1 form very stable nanostructures that match the dimensions found in the TEM image, we aimed to find distinctions by studying their infrared absorption in the amide I region. The absorption frequency and width of the amide bands is well-known to depend strongly on the vibrational coupling between adjacent moieties[49,50] and can be relatively easily measured experimentally. The experimental spectrum of monomer 1 (a single peptide strand with aromatic dithiol headgroup) is displayed in Figure A and exhibits a broad absorption centered on 1642 cm–1 as expected for a well-solvated, unstructured peptide. An additional contribution to the spectrum is present at 1594 cm–1, attributed to absorption of the C-terminal carboxylate group, which is deprotonated at pH 8. In comparison, the experimental spectrum of the stacked macrocycles displays a much narrower amide absorption red-shifted to 1622 cm–1 (Figure B), while the COO– absorption has lost some intensity and has shifted to 1590 cm–1. These changes in the spectrum are typical for short peptides self-assembling into β-sheet-like nanostructures as the amide moieties enter a much more well-defined hydrogen bonding network[51,52] and the pKa of the acid group is dramatically different in the assembled state, leading to protonation[53] or the formation of strong salt bridges. In all spectra an absorption band is present at ∼1673 cm–1, which is attributed to residual trifluoroacetic acid from the peptide synthesis. For clarity this band has been subtracted to produce the spectra in Figure (raw data in Figure S8).

Figure 4

Comparison of spectroscopic and theoretical results. All theoretical frequencies are red-shifted by 10 cm–1. Two different conformations (see Figure ) were modeled. (a) IR absorption of nonassembled peptide 1 in water (AmideIMaps: explicit SPC model, DFT: implicit continuum model). (b) IR absorption of fibers of 1. (c) UV–vis absorption of single fibers of 1. The calculated individual transitions are convoluted with a Gaussian curve with σ = 0.1 eV. Representative images of the involved transitions in 4 different bands are included as insets. (d) Temperature-dependent CD spectra of 1.The fit of the spectrum at T = 20 °C to Gaussian curves is indicated by the dashed line. The signal extrapolated to infinite temperature by a sigmoidal fit to the ellipticity as a function of temperature is shown in purple dots. (e) Molar ellipticity as a function of temperature for the 4 negative bands in the spectrum. The average temperature at which 50% of the CD signal at that energy has disappeared is indicated. (f) Experimental and calculated supramolecular CD signal for fibers of 1.

To test whether the proposed conformations match the experimental IR spectra, we used the MD trajectories to calculate the amide I absorption using a semiempirical vibrational exciton model combining frequency mapping with transition dipole couplings.[54]Figure A compares the calculated spectrum of monomer 1 on both the frequency mapping and the DFT level of theory with the experimental spectrum. The results of both methods show excellent agreement in terms of position and line width of the amide I absorption. The DFT result and experimental spectrum were used to calibrate the vibrational exciton result and indicated a horizontal redshift of 10 cm–1 needs to be applied to the frequency mapping results, which is in accordance with the average shift reported for the frequency mapping method.[54] A distinct absorption at 1608 cm–1 is present in both DFT and frequency map results as a consequence of the proximity of the first amide group to the aromatic ring, which lowers the absorption frequency. The corresponding area in the linear experimental spectrum is unfortunately too convoluted to conclude whether this peak is present (see Figure S8). Note that the COO– band from the DFT calculation was omitted from the spectrum for clarity (see Methods).

The vibrational exciton calculation on the fibers of 1 accurately predicts the narrowing and redshifting of the absorption (Figure B), indicating the proposed β-sheet packing is realistic. The overlap of the experimental spectrum is best with the cartwheel conformation, although the main absorption of the pairwise conformation is shifted only by 2 cm–1 and its width is only 3 cm–1 larger, both of which are within the error of the method. It is interesting to note that a strong but broad band is still present centered on 1645 cm–1 (see Figure S8), even in the spectrum of the fully converted library, which is also captured to some extent by the calculations. This indicates that structural heterogeneities, i.e., disordered or solvent-exposed amide groups are present in both the experiment and the calculations. The fact that these non-β-sheet vibrational modes are also captured by the model suggests that these heterogeneities are to some extent adequately represented in the MD simulations, although the limited size of the simulation (16 molecules ≈ 8 nm) could contribute and major defects in the fiber structure may be underrepresented. Interestingly, the spectrum of the double fibers displays a greater sheet-to-coil ratio (Figure S8), indicating less heterogeneities and disorder, which shows that slow oxidation leads to more well-ordered β-sheets.

As the calculations can be easily extended to a two-dimensional IR setup, performing this experiment could give more clarity in distinguishing the conformations in the future. Note that accurate DFT frequency calculations on systems this large (∼40 000 electrons) are currently inaccessible to most computational resources.

To further investigate the supramolecular structure within the fiber, we have calculated its electronic structure using the exciton method described above and in ref (34) to compare the UV–vis and supramolecular CD spectra of the two conformations to the experimental result. First, the semiempirically calculated spectra were calibrated using more accurate TD-DFT (CAM-B3LYP/cc-cpVDZ) calculations on 1 and its components. A good overlap was found between the spectra calculated at both levels of theory and the experiment (see Supporting Information, Figure S9), validating the use of CIS/ZINDO for the whole macrocycle system, for which the TD-DFT calculations are too computationally demanding.

Figure C shows a comparison of the experimental UV–vis absorption spectrum and the exciton calculation of the macrocycle fibers. Although the quantitative agreement between the calculated and experimental spectrum is somewhat disappointing, it is clear that the band structure of the transition is accurately predicted. Analysis of the location of the transition densities thus allows the determination of the origin of the various contributions to the spectrum; bands at 295, 238, and 210 nm originate from the ground-to-excited states transitions equivalent to the transitions to B2u, B1u and E1u states in benzene, respectively, but located on the aromatic core (see inset in Figure C). The first two of these transitions are no longer strictly forbidden due to symmetry breaking by the amide and disulfide groups. The phenylalanine ring in the peptide side chain only significantly contributes to the absorptions at 210 nm and below. Interestingly, it was found that the bands related to transitions on the aromatic headgroup are sensitive to the torsional angle of the amide group in both magnitude and position. Maximum intensity of bands at 310. 248 and 213 nm is achieved with a CCCN torsional angle around 0 or 180° (Figure S10). These bands are clearly present in the experimental spectrum, which strongly suggests the first amide bond of 1 is in the same plane as the aromatic ring, prohibiting, to some extent, H-bonding of the first amide group along the fiber axis. However, the aromatic rings are also tilted with respect to the fiber cross section, enabling some H-bond interactions as can be observed in Figure E,F. In any case, the MD results suggest that the flexibility and absence of a side chain in the first amino acid (glycine) are crucial in enabling H-bonding along the fiber axis by allowing the rotation of the following amide groups. This is confirmed by the notably different and unusual location (compared to common values for peptide secondary structure) of the glycine [φ, ψ] angles in the Ramachandran plot in comparison with the other amino acids in the peptide (Figure S11). Interestingly, the headgroup to glycine dihedral distribution of the cartwheel conformation in Figure S4 displays a clear asymmetry around 0°. This is a mark of the helicity of the fiber as this distribution is mirrored for d-amino acids and is clearly diminished in the pairwise conformation, indicating a less chiral structure as can be conveyed from Figure D,F. This tentatively suggests the cartwheel conformation as more realistic: from the TEM images in Figure a clear helical arrangement of the hexamer fibers was observed, suggesting a supramolecular chirality. Note that the UV–vis spectrum for the double fibers was qualitatively identical (Figure S12).

Previous work has demonstrated the sensitivity of the electronic CD signal to peptide sequence and macrocycle arrangement,[13] often exhibiting a typical β-sheet signal with a maximum around 194 nm and a minimum around 214 nm. This absorption shape was not clearly observed for the stacks of 1 as can be seen in Figure D, and as such, while the calculated UV–vis spectrum proved fairly insensitive to the different conformations studied, the chiral environment of the macrocycles is considered as a more promising marker.

The supramolecular component of the CD signal was extracted from temperature-dependent CD experiments by extrapolating the change in ellipticity upon heat-induced disassembly of the stacks (see Figure D–F and Methods). Considering the bisignated nature of CD transitions, the band positions in the supramolecular CD match the excitation energies observed in the UV–vis spectra of the disulfide headgroup and phenylalanine (Figure S9), which suggests the core of the macrocycle fibers is stacked in a chiral manner, not just on the level of the peptide strands, as seen for some other aromatic peptide amphiphiles.[32] The temperature at which half of the CD signal has disappeared (T1/2, Figure E) is 47 ± 2 °C for all bands in the spectrum, which is in agreement with the umbrella sampling MD results, where the binding free energy for peptide and aromatic components were found to be similar. Interestingly, we note that for the double fibers, even though the room temperature spectrum is qualitatively identical, T1/2 is notably higher (53 ± 1 °C, Figure S13), indicating that the fibers can be stabilized by bundling without major structural rearrangements, in line with the lower disordered component in the IR spectrum.

In addition to the experimental CD spectrum, the supramolecular CD signal was calculated using the excitonic model. Figure F shows the overlay of the calculated (averaged over 10 frames of the MD trajectory) and experimental supramolecular CD spectrum. As expected for a ZINDO calculation in vacuum, the calculated spectra are blue-shifted with respect to the experimental one. To facilitate their comparison, a shift of −0.3 eV was applied to the calculated spectra. The cartwheel and experimental spectra reasonably match, except at the high wavelength range where the transitions located on the aromatic core are poorly described. In this energy range, small changes in the conformation of the aromatic core, amide and disulfide groups can strongly affect the localization of the electronic transitions, and thus the CD spectrum. For the pairwise conformation, its small helicity leads to a rather flat CD signal. Furthermore, as the calculations were limited to relatively short stacks of 8 molecules (∼4 nm), it is possible that the calculations do not consider the right length scale as an origin of the chirality: from TEM images the period of the helicity of the double fibers was found to be 56 nm. As the CD signal is extremely sensitive to local disorder and heterogeneities,[35] undersampling may limit the applicability of the method on larger structures at this point. In any case, the fact that the CD signal does not match the archetypal β-sheet signal emphasizes again that care has to be taken when extrapolating protein database results to short peptides, especially when other chromophores are present.

The theoretical methods presented here have not been applied before to structures that have been self-assembled from randomly dissolved molecules in solution. The currently achievable time scale of MD simulations does usually not allow well-ordered structures, such as the structure in Figure , to form spontaneously. It is therefore possible that the results in this study miss a certain disordered component (i.e., defects) in the nanostructures, as suggested by the 1645 cm–1 band in the IR spectrum of 1. While this may cause additional broadening in the calculated spectra, the results are expected to be qualitatively unaffected as a consequence of the local nature of the excitonic coupling (both vibrational[55] and electronic[34]).

Conclusions

The nanoscale conformation of self-replicating macrocycles has been studied using electron microscopy, spectroscopy, molecular dynamics and semiempirical quantum chemical calculations. Single and double fibers can be prepared using fast or slow oxidation of dithiols to disulfide macrocycles respectively, where bundled fibers are more stable. The MD, IR, UV–vis and CD spectroscopy results support π–β-type interactions as the driving force behind the one-dimensional self-assembly, with approximately equal contributions from the peptide and the aromatic components to the binding free energy in the system studied here. The MD results suggest that extensive hydrophobic interactions between nonpolar residues in adjacent peptide chains occur. Furthermore, the simulations in which a ring was pulled off a stack suggests a mechanism for the reverse process of assembly (that was too rare to be captured in unbiased MD simulations) in which binding involves distinct steps in which the individual peptide chains of the incoming ring align with those in the stack and the aromatic core of the ring binds to the center of the stack. The multistep and cooperative nature of the process explains why stack nucleation is slow while elongation is faster. It also explains why stacks are of a uniform ring size, as peptide and headgroup alignment requires the number of peptide arms in the ring to match those in the stack, thereby favoring self-replication.

MD simulations also strongly suggest that, in contrast to the schematic representation in Figure and similar schemes in previous publications,[9,11,13] the aromatic core of the fibrils is collapsed and no water molecules are located within the core of the structure. The comparison of experimental data to the IR, UV and CD calculations on two hypothetical supramolecular architectures obtained from the MD simulations suggests that the “cartwheel” conformation, with all β-sheets’ hydrophobic faces in the same direction around the fiber axis, is the most likely supramolecular architecture, although more advanced experiments such as 2D-IR spectroscopy or solid-state NMR would be necessary to confirm this.

The semiempirical exciton methods to calculate IR, UV–vis and CD spectra from the result of the MD simulations have been applied to short peptides and their results are shown to be in good agreement with experimental spectra. These spectroscopic techniques are routinely used in the study of bionanostructures, but opening up the ability to link them directly with specific nanoscale architectures greatly enhances their usefulness in the field. The high level of structural insights achieved here by combining experimental and theoretical approaches demonstrates their usefulness to investigate other bioinspired self-assembling molecules.

Methods

Molecular Dynamics Simulations

All MD simulations were performed using the GROMOS 54a8 force field,[56,57] with in-house developed parameters for the aromatic disulfide headgroup (see Supporting Information for parametrization strategy). Simulations were started from a straight (achiral) stacked structure. SPC water[58] was added in a layer of at least 1 nm thickness around the fiber and enough chloride ions were added to neutralize the system. All explicit hydrogen atoms were converted to virtual sites and all molecular bonds were constrained in production runs using the LINCS algorithm,[59] except for water molecules which were constrained by the efficient SETTLE algorithm.[60] Additionally, the center-of-mass motion of the solute (both translation and rotation) were removed every 100 time steps. After a short minimization of 5000 steepest descent steps, all runs were performed with time steps of 2.0 fs in the NPT ensemble using the velocity-rescaling thermostat[61] to separately couple the peptide macrocycles and the rest of the system (water and ions) to an external bath at 298 K with a coupling constant of 1.0 ps. The pressure was kept at 1.0 bar using the Berendsen barostat[62] with a coupling constant of 1.5 ps. A Barker–Watts reaction field with εRF = 62 was used to treat long-range electrostatic interactions with Coulombic and van der Waals forces cut off at 1.4 nm. All MD simulations were run using GROMACS 4.6.7.[63] Trajectories of 100 ns were generated in two steps: during the first 50 ns distance restraints with a force constant of 1000 kJ/mol were applied to keep Cα atoms of nearest neighbors in the other peptide strands at a constant distance of 0.48 nm. These restraints were then removed for the second period of 50 ns.

Umbrella sampling simulations to determine the binding free energy of truncated macrocycles were performed on a trajectory created by pulling the top molecule by its center of mass (COM) from a pre-equilibrated fibril of either 2 or 8 macrocycles along the fibril axis, while keeping position restraints (1000 kJ/mol) on all Cα atoms of the second macrocycle in the stack. The pulling trajectory was divided into windows with a 0.12 nm COM spacing and an NPT simulation was performed for 10 ns within each window. The free energy profile of the dissociation process was calculated using the weighted histogram analysis method (WHAM) as described in ref (46) with Bayesian bootstrapping for the error estimate. As truncated peptide dimer stacks were not stable in the simulation, their binding free energy was determined by counting the number of dimers (dimers are defined by having any two atoms within a distance of 0.4 nm) during a 1 μs NVT simulation of two peptides in a cubic 658 nm3 box of SPC water, pre-equilibrated in the NPT ensemble for 500 ns. The binding free energy was determined from the Boltzmann ratio of monomers to dimers as performed in ref (64).

Vibrational Exciton Analysis

Calculations of the vibrational frequency of all amide groups was performed on the result of a 1 ns MD simulation added to the end of the simulations described above, saving the system’s coordinates (including water molecules) every 20 fs, thus averaging over 50 000 frames. The vibrational Hamiltonian was created using the AmideIMaps program v1.0.0[54] with the Skinner frequency map[28] and Torii’s transition dipole coupling (TDC),[65] as recommended by Bondarenko and Jansen in combination with GROMOS force fields.[54] The NISE2A program was employed to calculate the linear absorption as determined by the Fourier transform of the first-order response function.[27,66,67] Ensemble averaged spectra were obtained by evaluating the response functions from 5000 starting configurations spaced by 200 fs along the trajectories. All spectra were red-shifted by 10 cm–1 to maximize spectral overlap, in accordance with the findings of Bondarenko and Jansen for the combination of GROMOS54a7, Skinner map and TDC. Note that no posthoc broadening was applied to the results; the width of the calculated spectra is the natural line width originating from the vibrational motion of the molecules as sampled by the MD simulation.

DFT Vibrational Analysis

DFT calculations on the zwitterionic monomer (lysines protonated, C-terminus deprotonated) were performed using Gaussian09 rev. D.01[68] employing the Grimme’s dispersion-corrected B97-D functional[69] with def2-SVP basis set[70] to optimize the geometry (using loose convergence criteria) and calculate the harmonic frequencies. An empirical frequency scaling factor of 0.965 was applied to correct for basis set truncation and anharmonic effects and the calculated absorption frequencies were convoluted with a 12 cm–1 Gaussian line shape function to simulate the effect of the dynamic solvent environment. Vibrations located the COO– functionality were omitted from the result as their frequency is not accurately predicted using this method (calculated: 1611 cm–1, experimental: 1593 cm–1) and obscures the result, see also ref (51).

UV–Vis and CD Calculations

The geometric structures generated from the simulation runs above (10 snapshots from the last 5 ns of the simulation for cartwheel and pairwise conformations) were converted to an explicit hydrogen representation using the default H–C distances from the CHARMM27 force field.[71,72] These full coordinates were used as input for quantum-chemical calculations at the ZINDO/CIS level of theory using the Gaussian09 program. An active space of 200 occupied and unoccupied orbitals and the first 500 excitations were taken into account (1000 for a single macrocycle). The calculations were performed on stacks of 8 molecules (instead of 16) due to computational memory limitations. The results of these calculations served in turn as the input for a phenomenological Frenkel exciton model to predict the spectroscopic properties of the helical supramolecular structures as described in detail previously.[34,35,73] Absorption energies were convoluted with a 0.10 eV wide Gaussian function. Note that signals <200 nm are less reliable due to the cutoff at 500 excitations per molecule. For calibration and comparison, TD-DFT spectra of the monomer (in the carboxylic acid form) were calculated using Gaussian09 employing the CAM-B3LYP functional with a coupled-cluster polarized valence double-ζ basis set. All UV–vis and CD calculations were performed in vacuum. A redshift of 0.3 eV was necessary to correct for the absence of solvent and method artifacts as confirmed by calculations and experiments on various model compounds of 1.

Library Preparation and Analysis

Peptide 1 was purchased from Cambridge Peptides Ltd. and used as received (purity >90%). The libraries were prepared by dissolving the compound 1 (10 mM) in borate buffer (50 mM, pH 8.2). For preparation of single fibers only, the library was partially oxidized by the addition of 0.8 equiv of NaBO3 (Sigma-Aldrich, from a 20 mM aqueous solution). All libraries were left agitated in the presence of oxygen from the air until the hexamer was the main species (>95%) as confirmed by UPLC-MS. Immediately after preparation, the pH of the solution was adjusted by the addition of 1 M KOH solution such that the final pH was 8.2. All library experiments were performed at ambient temperature. A small aliquot of each sample was removed to another vial and diluted 20 times with doubly distilled water prior to HPLC or UPLC(-MS) analysis.

The buffer was prepared from anhydrous borax (Fluka) and boric acid (Merck Chemicals) dissolved in doubly distilled water from in-house double distillation facilities. Libraries were prepared in clear HPLC glass vials (12 × 32 mm) closed with Teflon-lined snap caps purchased from Jaytee. Library solutions were stirred using Teflon coated microstirrer bars (2 × 2 × 5 mm) obtained from VWR. Samples were stirred on a Heidolph MR Hi-Mix D magnetic stirrer at 1200 rpm.

Fourier Transform Infrared Spectroscopy

Libraries of 1 were freeze-dried and resuspended in D2O (99.9% D, Sigma) to the same concentration (1.7 mM in hexamer) at pD 8.0. The sample was kept overnight to mature and was then sandwiched between two CaF2 plates separated by a 50 μm PTFE spacer. Monomer samples were prepared by dissolving peptides in 0.050 M deuterated borate buffer at pD 8.0 and measuring within 5 min of preparation to prevent oxidation of the thiol groups. FTIR spectra were obtained from an average of 32 spectra on a PerkinElmer Spectrum Series 400 spectrophotometer with a resolution of 4 cm–1. The resulting spectra were corrected for D2O and water vapor absorption.

UV–Vis Spectroscopy

UV–visible measurements were performed using a Jasco V-650 UV–visible spectrophotometer at 298 K. A quartz cuvette with 1 cm path length was used for the measurements. The absorption spectra were recorded from 190 to 850 nm at a 200 nm/min scanning rate. Monomer and hexamer samples were measured at a concentration of 1.0 × 10–5 M in pentapeptide. It was confirmed the spectra did not change qualitatively at higher concentrations.

Temperature-Controlled Circular Dichroism Spectroscopy

Spectra were obtained on a Jasco J-810 spectrometer with a Peltier temperature controller. Samples were diluted to 5.4 × 10–5 M in pentapeptide were heated/cooled from 283 to 343 K in steps of 1 K (5 K for single fibers) at a rate of 0.1 K/min and kept for 5 min at every temperature before measuring. Spectra were obtained as averages of 3 measurements from 190 to 300 nm with a scanning speed of 200 nm/min and a bandwidth of 1 nm. The spectra as a function of excitation energy were fitted to Gaussian bands around 6.5, 6.2, 5.8, 5.2, and 4.6 eV using the global fit routine in Igor Pro with free parameters for the band position, width and amplitude. Ellipticities as a function of temperature thus obtained were fitted to a sigmoidal curve to extrapolate the CD signal lost at high temperature, which equals the supramolecular component of the CD signal. Data points for which the chi-square of the fit was larger than 100 were considered not converged and removed from the graph (region λ < 198 nm).

Transmission Electron Microscopy

Samples were diluted to approximately 0.80 mM using doubly distilled water. A small drop (∼5 μL) of sample was deposited on a 400 mesh carbon-coated copper grid (Agar Scientific) and blotted on filter paper after 30 s. The sample was stained twice using 4 μL of 2% uranyl acetate solution. The grids were imaged in a Philips CM12 electron microscope operating at 120 kV using a slow scan CCD camera. Fiber widths were measured at ≥200 locations across 5 images (7 for double fibers), from the inner layer of the dark contrast layers observed.