Ion Effect and Metal-Coordinated CrossLinking for Multiscale Design of Nereis Jaw Inspired Mechanomutable Materials Chia-Ching Chou,† Francisco J. Martin-Martinez,† Zhao Qin,† Patrick B. Dennis,‡ Maneesh K. Gupta,‡ Rajesh R. Naik,‡ and Markus J. Buehler*,† †
Laboratory for Atomistic and Molecular Mechanics, Department of Civil and Environmental Engineering, Massachusetts Institute of Technology, 77 Massachusetts Avenue, Cambridge, Massachusetts 02139, United States ‡ Air Force Research Laboratory, Materials and Manufacturing Directorate, Wright-Patterson AFB, Ohio 45433, United States ABSTRACT: The Nvjp-1 protein is a key component in the jaws of Nereis virens, a species of marine worm. It contains over 25 mol % of histidine, which is believed to play a key role in the metal-coordinated cross-linking responsible for the structural stability and exceptional mechanical performance of the worm jaw. Understanding the nanoscale mechanism behind this cross-linking and its pathway in affecting the macroscopic mechanical behavior of the material is crucial to develop bioinspired mechanomutable materials based on Nvjp-1. Here, we use a combination of multiscale modeling and experimental synthesis to understand the behavior of this heterologous-expressed protein from the nano- to the macroscale. We have built a bottomup molecular-based model, which includes electronic-based density functional theory calculations, atomistic simulation of the nanoscale properties with replica exchange molecular dynamics, and an elastic network model for describing the macroscale behavior at different pHs. This multiscale modeling supports the experimental synthesis of a photo-cross-linked Nvjp-1 hydrogel by proving both the nanoscale mechanisms and mechanical behavior predictions. Our theoretical results agree well with the experimental observations, showing that Nvjp-1 forms a more compact structure in the presence of Zn2+ ions with a suitable pH environment, leading to the formation of more stable intramolecular metal-coordinated cross-links. These metal-coordinated cross-links induce nanoscale aggregation of Nvjp-1, which is responsible for the hydrogel contraction observed in experiments and predicted by the model. KEYWORDS: Nvjp-1 protein, metal coordination, Nereis jaw, molecular dynamics, network model
T
electrons. Even though the coordination spheres in the metal atom have a well-established structure, and some of them create bonds almost as strong as covalent ones, they have been shown to serve as a flexible connection between ligand and ion.8,9 Furthermore, the strength of the bond usually decreases in lowpH conditions and increases with pH. Most importantly, unlike most covalent bonds, whose rupture is usually irreversible, the dative nature of metal-coordinated bonds makes them able to re-form after rupture, which enables the ligand and ion exchange process. Nereis worms have taken advantage of this mutable strength of metal-coordinated bonds in different chemical microenvironments to develop an intriguing Nvjp-1 mutability in aggregation processes and hydrodynamics proper-
he Nvjp-1 protein is a key component of Nereis jaw, which is composed of roughly 90% protein, 8% halogen atoms, and 1% Zn2+ ions.1−3 Surprisingly, although most of the jaw of this species of worm is composed of organic matter, its hardness and stiffness are comparable to that of human dentin (∼1−2 GPa of hardness and 10−20 GPa of stiffness). More than 25 mol % of the Nvjp-1 amino acid sequence is histidine, which commonly forms metal-coordinated bonds with ions such as Cu2+ or Zn2+. In fact, previous studies have pointed out that this formation of metalcoordinated cross-links between histidine-rich Nvjp-1 and Zn2+ ions dominates the structural stability of the protein, and it is responsible for the excellent mechanical properties of the Nereis jaw.1−7 As it is well known, metal-coordinated bonds are dative chemical bonds that occur when one of the species, i.e., ligand, donates a lone pair of electrons to a metal ion, usually a transition metal, with empty orbitals to locate such © 2017 American Chemical Society
Received: November 23, 2016 Accepted: February 6, 2017 Published: February 6, 2017 1858
DOI: 10.1021/acsnano.6b07878 ACS Nano 2017, 11, 1858−1868
Article
www.acsnano.org
Article
ACS Nano
Figure 1. Molecular models of histidine−metal complexes with Na+ and Zn2+. All the structures are fully optimized with B3LYP-gCP-D3/631G*. The energies of each molecule relative to the most stable structure (considered as zero) are shown in parentheses for the three sets of molecules. Bond distances (Å) between metal ions and histidine functional groups are indicated close to the dashed lines in the figure. (A) Systems composed by two histidine molecules and one Na+ ion, H2Na-I, H2Na-II, and H2Na-III, show three different conformations studied. The relative energy of H2Na-II and H2Na-III with respect to the most stable one (H2Na-I) is shown in parentheses. (B) Systems composed by two histidine molecules and one Zn2+ ion, H2Zn-I, H2Zn-II, and H2Zn-III, show the same three different conformations studied. The relative energy of H2Zn-I and H2Zn-III with respect to the most stable one (H2Zn-II) is shown in parentheses. (C) Systems composed by three histidine molecules and one Zn2+ ion, H3Zn-I and H3Zn-II, show two different conformations studied. The relative energy of H2Zn-II with respect to the most stable one (H2Zn-I) is shown in parentheses. A side view for H3Zn-I is also shown for clarity.
ties. As an example, as pH increases in a solution of purified Nvjp-1, the sedimentation profile from ultracentrifugation indicates multiple species of high-order oligomers with decreased sedimentation coefficients,6 which suggests a significant change in the molecular shape of the protein in response to the change in pH. Despite its key role in the Nereis jaw, relatively little is known about the mechanisms of the mutable metal-coordinated bonds occurring in Nvjp-1 proteins. Particularly, more experimental and theoretical research is needed to get further insights into
the different behavior of this protein in different environments. A predictive model that encompasses different length scales is critical to reach this understanding and to probe the role of Nvjp-1 nanostructure and metal coordination on the mechanical properties of the Nereis jaw. Such a multiscale model will also help to clarify the origin of the abovementioned mutability observed in the behavior of Nvjp-1 and to establish the link between the protein’s chemical structure, abundance of metal-coordinated cross-links, nanoaggregate formation, and large-scale properties. It will also benefit the 1859
DOI: 10.1021/acsnano.6b07878 ACS Nano 2017, 11, 1858−1868
Article
ACS Nano Table 1. Binding Energiesa
design of Nereis jaw-inspired materials in current studies8,10 and will help to exploit natural mechanisms for materials bioengineering and many applications. Moreover, developing a multiscale model of this protein, together with the subsequent predictive theoretical framework, enables further bottom-up studies on the mechanical properties of similar biological structures. In this work, we focus on the effects of Zn2+ ions on the structure and mutability of the Nvjp-1 protein at the nanoscale, and we investigate the role of the metal coordination in the cross-linking mechanisms by using density functional theory (DFT) calculations and replica exchange molecular dynamics (REMD) simulations. We also predict the large-scale composite material behavior by using large-scale coarse-grained models based on atomically trained elastic networks. We couple this with the experimental synthesis of a photo-cross-linked Nvjp-1 hydrogel. This combination of theory and experimentation provides fundamental insights, from a bottom-up perspective, into the mutable mechanical properties of the Nvjp-1 protein. This is one of several mechanisms by which the observed contraction can take place, being Manning-type charge condensation another mechanism that might be expected to dominate in some layer at the periphery of the gels. These results can be extrapolated to other similar proteins with metalcoordinated bonds and provide a theoretical basis for designing mechanomutable polymer materials from the most fundamental scale and up.
complex
binding energy [kcal/mol]
H2Na-I H2Na-II H2Na-III H2Zn-I H2Zn-II H2Zn-III H3Zn-I H3Zn-II
−50.88 −53.45 −49.40 −231.68 −231.45 −226.23 −247.89 −217.72
a
The values are calculated as the difference between the energy of the complex and the addition of the energies of the different components separately.
amount of nitrogen atoms involved in the interaction (H2ZnII) is also the most stable. This suggests that the interaction between histidine residues and Zn2+ largely occurs through nitrogen atoms in the protein. In fact, this is consistent with the primary structure of the protein, where oxygen atoms are located closer to the peptide bond and therefore in a more isolated position regarding interactions with ions in the environment. Further insight into the fundamental nature of these interactions was acquired by analyzing the frontier molecular orbitals of these complexes. Figure 2 shows the highest occupied molecular orbital (HOMO) and the lowest unoccupied molecular orbital (LUMO) for all the histidine− metal complexes with Na+ and Zn2+. Figure 2A and B show that the HOMO and the LUMO are mostly localized in the histidine molecules rather than in the ions. Actually, the shape and location of these orbitals do not change with respect to the isolated histidine molecule (Figure 2C). This suggests a mostly electrostatic interaction between the positively charged ions and the lone pairs of electrons of either oxygen or nitrogen atoms. On a separate note, it is worth mentioning that those complexes that were shown to be less stable (H2Na-III and H2Na-III) are also the ones with a symmetry break in their molecular orbital shape. In these cases, the HOMO and the LUMO orbitals are localized only in one of the histidine molecules, and this symmetry break in the electronic structure destabilizes the complex. In the case of three-coordinated complexes, only Zn2+-based structures are shown in Figure 2D. In this case, the HOMO is degenerated, and both HOMO and HOMO−1 are shown. It is remarkable that in a three-coordinated arrangement the LUMOs are localized both at the ion and at the histidine molecules. This implies a stronger electronic interaction, in addition to the electrostatic one, where electrons are located on both sides in the molecular system. 2.2. Replica Exchange Molecular Dynamics (MD) Simulations on Nvjp-1 Protein. To understand how these metalcoordination interactions affect the protein’s folding and aggregation structures and thereafter their mechanical functions, we perform atomistic MD simulations and investigate how ions affect the molecular conformation. Four lowestenergy representative structures were obtained from REMD simulations in implicit solvent, as shown in Figure 3. According to our simulation results, the Nvjp-1 protein has a primary turn and random coil configurations and forms small clusters of beta-sheet and helix secondary structures. Different secondary structure contents are illustrated in Figure 3A to D. Both N(residue 1−174) and C-terminals (residue 175−381) contain
RESULTS AND DISCUSSION Density Functional Theory Calculations on CrossLinking Sites. The bottom level of this multiscale model implies a complete description of the electronic structure and the electron density at the cross-linking sites. To get a deep understanding of the fundamental mechanisms occurring at the nanoscale, we perform DFT calculations on model systems that represent the cross-linking points mediated by metal ions. To evaluate the effect of Zn2+ in comparison to Na+, we analyze different complexes of histidine with Zn2+ and Na+ in implicit water. Figure 1 shows the DFT fully optimized geometries for the histidine−metal complexes under study. Different arrangements between ions and histidine molecules have been explored to identify the most favorable orientations, as well as to identify preferable functional groups for the interaction. Both complexes with two and three histidine molecules have been explored to evaluate the most suitable cross-linking arrangement between protein terminals, i.e., two- or threecoordinated. In Figure 1, H2Na-I, H2Na-II, and H2Na-III are Na+ complexes, while H2Zn-I, H2Zn-II, and H2Zn-III are Zn2+ complexes. Zn2+ complexes always exhibit shorter bond distances. This implies stronger interaction, which agrees with the contraction observed in the experiments upon exposure to Zn2+, as discussed below. These experimental results also agree with the DFT-calculated binding energies that are shown in Table 1. Binding energy values clearly indicate that the formation of Zn2+ complexes further stabilizes the structure in comparison to Na+ complexes. In fact, the binding energies in Zn2+ complexes are more than 1 order of magnitude larger than those for Na+ complexes. These results support an ion exchange mechanism, where Na+ ions interacting with histidine terminals would be replaced by Zn2+ ions upon addition of Zn2+ to the NaCl solution. Additionally, DFT results also show that in the case of Zn2+ complexes the molecular conformation with the highest 1860
DOI: 10.1021/acsnano.6b07878 ACS Nano 2017, 11, 1858−1868
Article
ACS Nano
Figure 2. Molecular orbitals of histidine−metal complexes with Na+ and Zn2+ represented at 0.05 isosurface value. All the structures are fully optimized with B3LYP-gCP-D3/6-31G*. Highest occupied molecular orbital (HOMO) and lowest unoccupied molecular orbital (LUMO) are shown. (A) HOMO and LUMO for two-histidine-Na+ complexes, i.e., H2Na-I, H2Na-II, and H2Na-III. (B) HOMO and LUMO for twohistidine-Zn2+ complexes, i.e., H2Zn-I, H2Zn-II, and H2Zn-III. (C) HOMO and LUMO for the histidine molecule, for comparison with the complexed situation. (D) HOMO and LUMO for three-histidine-Zn2+ complexes, i.e., H3Zn-I and H3Zn-II.
from the synthesized protein hydrogel cannot utilize evaporation as do other structural protein polymers, like spider and silkworm silks. In fact, the initial driving force for removing bulk water from the hydrogel is likely due to charge screening of protonated histidine by the anion. The type of cation present in the salt seems to play a role in the degree of contraction experienced by the hydrogel, where zinc treatment consistently yields a more compact structure compared to that seen with sodium. We observed this contraction behavior upon addition of as little as 10 mM NaCl, and it reached nearly maximum contraction at a concentration of 200 mM NaCl (Figure 4B). In order to discern between contraction caused by monovalent and divalent cations, fully swollen Nvjp-1 bars were incubated in 200 mM NaCl, 100 mM ZnCl2, and 200 mM ZnCl2 (Figure 4C). While contraction was observed with all salts, the bars incubated in both 100 and 200 mM ZnCl2 contracted to a greater extent than those in 200 mM NaCl. Bars were incubated in two concentrations of ZnCl2 to account for the ionic equivalents that can drive osmotic drying of hydrogels. As a result, the additional 6% contraction observed between the 200 mM NaCl and the 100 mM ZnCl2 is most likely due to the coordinate cross-linking from the Zn2+ cations. In order to support and better understand these experimental results, we calculated the solvent-accessible surface area (SASA) from the REMD simulations at different ion concentrations
roughly the same number of clusters of beta-sheet and helical structures, although the C-terminal exhibits relatively longer clusters of helix. This difference might be due to the lower content of glycine amino acids in the C-terminal domain (30%) with respect to the N-terminal domain (42%). On the basis of these representative structures obtained from implicitly solved REMD, we perform explicit-solvent REMD simulations in order to understand the mechanism of mutability of the Nvjp-1. These explicit-solvent REMD simulations with various ion concentrations probe the nanostructure of the Nvjp-1, as well as its mutable properties at the nanoscale. The mutability of the Nvjp-1 protein was already observed in previous studies, where significant Nvjp-1 aggregation occurred at low Zn2+/protein molar ratios, with aggregation slightly increasing at a Zn2+/protein molar ratio of 50.6 Drop-cast films made from the purified Nvjp-1 were not stable in neutral salt solutions, breaking up into large protein aggregates. Synthesis of Photo-Cross-Linked Hydrogels in Experiments and Comparison with Simulations. Here, we show that photo-cross-linking of purified Nvjp-1 creates a stable hydrogel that can be made to significantly contract or expand upon exposure to different ions and concentrations (see Figure 4). Since Nereis virens is a marine organism, the dehydration of this nascent worm jaw material during biosynthesis must occur in an aqueous environment. Thus, the removal of bulk water 1861
DOI: 10.1021/acsnano.6b07878 ACS Nano 2017, 11, 1858−1868
Article
ACS Nano
Figure 3. Four representative structures from implicit solvent REMD simulation. The coloring is based on the secondary structure of each amino acid using the STRIDE algorithm in VMD.15
Figure 4. Ionic effect on Nvjp-1 protein observed in experiments. (A) A photo-cross-linked Nvjp-1 hydrogel was swollen at low pH (−NaCl) before being exposed to 500 mM NaCl in water (+NaCl). The scale bar is 1 cm. (B) Contraction behavior of photo-cross-linked Nvjp-1 hydrogels after exposure to increasing NaCl concentrations. The degree of contraction of a hydrogel bar (A/A0) is expressed as a fraction of the final contracted dimensions (A = length × width) and its dimensions after placement in water without salt (A0). (C) Contraction of photocross-linked Nvjp-1 hydrogels with NaCl and ZnCl2 salts. The degree of hydrogel contraction is expressed as in (B), where the final contracted dimensions of a hydrogel bar are expressed as a fraction of the same bar swollen in a 10 mM HCl solution. Representative hydrogel bars are shown above in the fully contracted state where the surrounding dotted box represents the original swollen dimensions of the hydrogel.
(Zn2+/protein ratio = 0, 20, 50). The results are shown in Figure 5A. The SASA of the Nvjp-1 protein decreases as the Zn2+/protein ratio increases, which agrees with the Nvjp-1 solubility measured in experiments.6 From the Zn2+/protein ratio from 0 to 20, the SASA value decreases from 28 000 Å2 to approximately 26 250 Å2, indicating less solvent exposure of structure-buried residues. While the Zn2+/protein ratio increases from 20 to 50, the SASA decreases only slightly at the Zn2+/protein ratio of 50. We further calculate the radius of gyration to measure the size of the Nvjp-1 in these different conditions. Figure 5B shows the radius of gyration of Nvjp-1 as a function of Zn2+/protein ratio. The radius of gyration of the Nvjp-1 protein decreases as the concentration of Zn2+ ions increases. These results suggest that more compact structures are formed in the presence of Zn2+ ions, which is indeed consistent with the SASA results.
Looking deeper into the ion effect on metal-coordinated cross-links of the Nvjp-1 protein, Figure 6 shows that the number of histidine residues and Zn2+ ions involved in metal coordination increases as the Zn2+/protein ratio increases. However, the formation of metal-coordinated cross-links is almost saturated at lower Zn2+/protein ratios, and both the number of histidine residues and the number of Zn2+ ions rise only slightly at a Zn2+/protein ratio of 50. This pattern implies that more intramolecular cross-links are formed as the Zn2+/ protein ratio increases and also that the Zn2+ ions stay close to the center of histidine residues to form stable cross-links in a chain. However, the number of cross-links does not look to increase much more when higher Zn2+/protein ratios are reached, because saturation is achieved. In addition to the SASA and radius of gyration calculations already discussed, these results on the formation of a compact structure with intrachain 1862
DOI: 10.1021/acsnano.6b07878 ACS Nano 2017, 11, 1858−1868
Article
ACS Nano
Figure 5. Solvent-accessible surface area (SASA) and radius of gyration of the Nvjp-1 as a function of Zn2+/protein ratio. (A) The SASA of the Nvjp-1 protein decreases as the Zn2+/protein ratio increases, showing the same trend in the solubility of the Nvjp-1 protein measured in experiments. (B) The radius of gyration of the Nvjp-1 protein decreases in the presence of Zn2+ ions.
Figure 6. Number of (A) histidine amino acids and (B) zinc ions forming in metal coordinate complexes of Nvjp-1 protein at varied ion concentrations.
Figure 7. Simulation snapshot of intramolecular cross-link forming in a metal-coordinate complex. Zn2+ ions and histidine amino acids are highlighted, and the protein is colored based on the secondary structure. Water molecules are not shown.
cross-links provide a plausible explanation for the contraction behavior observed in the experiment. Note that because we use
a single protein in the simulation, only the intramolecular metal-coordinated cross-links are described. 1863
DOI: 10.1021/acsnano.6b07878 ACS Nano 2017, 11, 1858−1868
Article
ACS Nano
Figure 8. Elastic network model of Nvjp-1 material and its mechanical response in different pH conditions. (A) The computational model is a bead−spring model with an fcc lattice. The equilibrated geometry and stiffness of the material is reflected by the spring constant (2k) and the equilibrium length constant of the lattice structure (a), and both of their values are given as a function of the pH value in simulations. (B) Simulation results of a material ribbon under two different pH conditions, in comparison with the experimental observations of the Nvjp-1 hydrogels. One Nvjp-1 hydrogel was swollen in 10 mM HCl, while the other was contracted in 20 mM EDTA, pH 8.
Figure 9. Computational prediction of the mechanical behavior of Nvjp-1 material with top and bottom surfaces been subjected to a pH gradient. (A) Given by the snapshots taken during the molecular dynamics simulation, a long narrow ribbon with top surface subjected to high pH and bottom surface subjected to low pH forms the coiled structure with its appearance in analogy to a coil in plants (as inserted). (B) For the same simulation with different initial geometry, a round disk made of Nvjp-1 material with top surface subjected to high pH and bottom surface subjected to low pH folds and forms a cypraea-like geometry (in analogy, as inserted). It is shown by our simulations that Nvjp-1 can be used to design active structures that can be used to sense a pH gradient.
Analysis of Protein Sequence (SAPS)11 as repetitive motifs in the Nvjp-1 protein, while some repetitive histidine-rich peptides have been characterized to form metal-coordinated bonds in proteins.12−16 Since glycine amino acids add local structural flexibility to the protein, the histidine amino acids close to them are allowed to form cross-links without causing
On the basis of the formation of intrachain cross-links in the Nvjp-1 protein, we further analyze the relationship between the Nvjp-1 sequence and the location of metal-coordinated crosslinks. The analysis indicates that the cross-links mainly occur in the glycine- and histidine-rich segments, namely, GHHG and GHGH. These segments have been identified by Statistical 1864
DOI: 10.1021/acsnano.6b07878 ACS Nano 2017, 11, 1858−1868
Article
ACS Nano
screening of protonated histidine by the anion. The type of cation present in the salt seems to play a role in the degree of contraction experienced by the hydrogel, where Zn2+ treatment consistently yields a more compact structure compared to that seen with Na+. Our DFT calculations show a stronger interaction at cross-linking sites arising from the presence of Zn2+ ions, for which a three-coordinated cross-linking is slightly preferred, providing the most fundamental insight of the mutable mechanics of such biopolymer material. In addition, our atomistic model of the Nvjp-1 protein provides a better understanding of this phenomenon and predicts that the Nvjp1 has a different mechanical response in different solutions with various ion concentrations. Structurally and mechanically tunable properties of Nvjp-1 are demonstrated with this model and compared against the experiments. We have calculated SASA of Nvjp-1 with the atomistic model and have found that the SASA of Nvjp-1 decreases as the Zn2+/protein ratio increases. We have also calculated the radius of gyration to measure the size of Nvjp-1, and the results show that a more compact structure is formed in the presence of Zn2+ ions. These results show the consistent trend for solubility as with the values measured in the experiments. Our coarse-grained model of Nvjp-1 provides mesoscopic insights of its mechanical response, which demonstrates that it can adapt to shapes very different from its initial states by generating large mechanical deformations, providing opportunities to predict and design composite materials of active mechanical functions. Our results suggest that metal-coordinated cross-links play a significant role in achieving the characteristic properties of this protein material. Nvjp-1 provides an excellent example of mechanomutable protein material to investigate. By integrating multiscale analysis and experimental validation, we can learn the fundamental material science concepts at different scales from the quantum level of bond formation to the atomistic scale protein folding and mesoscopic scale swelling and contraction in different chemical conditions. Our multiscale modeling provides the most comprehensive understanding of the formation, mechanical behavior, and the tunable function of such protein materials under changing conditions. Inspired by Nvjp-1, we will use this computational modeling tool for rational polymer material design with desired tunable mechanical functions, which provides a rational way to fully program the mechanical response of the material in different conditions, which is missing in most previous studies. This study provides fundamental insight into mutable mechanical properties of Nvjp-1 and other proteins with metal-coordinated bonds from a bottom-up perspective. It also provides an opportunity to synthesize biomaterials with mutable properties in the presence of metal-coordinated cross-links. Additionally, these pH-sensitive materials can be used for the design of bioinspired adaptive materials for chemical sensing in cyber-physical systems or for active control of the deformation and motion of soft robotics, as it is shown in our large-scale simulations from which plant coil behavior and cypraea-like geometries are generated.
noticeable steric effects. A similar feature of glycine- and histidine-rich peptides binding to metal ions has been already described.10,15 Figure 7 shows a comparison of our simulation with the proposed model of cross-links in the Nereis jaws. The simulation snapshot clearly shows that Zn2+-mediated crosslinks connect segments in the histidine-rich protein, which implies that the metal-coordinated bond is indeed responsible for stabilizing the material and further enhancing the mechanical properties of the Nereis jaws. Elastic Network Model of Nvjp-1 Matching Experimental Behavior. Finally, we have developed an elastic network model of Nvjp-1 hydrogels to simulate their mechanical response in different pH conditions. Figure 8 shows the definition of the model and its performance against the experimental results. The simulation results of a material ribbon under different pH conditions are shown in comparison with the experimental observations of the Nvjp-1 hydrogels. Our model is able to simulate contraction occurring in experimental tests for two different pH values (pH 2 and 8), and these pH effects can be used for designing foldable and adaptive materials. Additionally, to complete the multiscale model and to provide some potential applications arising from these materials, we study how their functionality depends on the shape, structure design, and environmental conditions. Figure 9 shows an MD simulation of the behavior of a pH-sensitive Nvjp-1 material under different environmental conditions. Figure 9A shows the rolling of a long narrow ribbon with the top surface subjected to high pH and the bottom surface subjected to low pH. The coiled structure arising from this design resembles the coil occurring in certain plants. On the other hand, a similar material with a different initial geometry gives completely different shapes. In this case, Figure 9B shows a round disk made of Nvjp-1 with top surface subjected to high pH and bottom surface subjected to low pH, which folds into a cypraea-like geometry. It is important to note that while the natural shapes mimicked are obtained by differential growth in nature, those in the simulations are obtained by differential contraction. These two examples emphasize the potential of our bioinspired design and Nvjp-1 hydrogels to design active structures that can be used for pH gradient sensing, pH dynamic materials, soft robotics, and many other applications arising from these functionalities. Also, this model enables the development of a predictive framework for further design of Nvjp-1-based materials in the future.
CONCLUSION We provide a complete multiscale modeling framework and experimental synthesis of bioinspired hydrogels based on Nvjp1 protein. The jaws of N. virens are condensed proteinaceous structures whose mechanical properties are tuned by coordination with transition metal cations, particularly Zn2+. N. virens is a marine organism, and the dehydration during biosynthesis must occur in an aqueous environment. Since the removal of bulk water from the synthesized protein hydrogel cannot utilize evaporation, it is possible that the contraction behavior observed in this study mimics an as yet unknown natural mechanism for dehydration of worm jaw material as it matures into a final functional form. The contraction of acidswollen Nvjp-1 hydrogels is not dependent on the type of cation present in the salt. In fact, the initial driving force for removing bulk water from the hydrogel is likely due to charge
METHODS Replica Exchange Molecular Dynamics Simulation. REMD simulations17 are used to investigate the formation of protein structures in various environments. REMD is an effective tool that allows us to identify a native protein structure by overcoming kinetic 1865
DOI: 10.1021/acsnano.6b07878 ACS Nano 2017, 11, 1858−1868
Article
ACS Nano
Figure 10. Illustration of simulation protocol for identifying the structure of Nvjp-1 in this study. trapping in local energy minima in the simulation. The simulation approach is illustrated in Figure 10. Primary protein structures are created with the amino acid sequence of the Nvjp-1 protein.6 To simulate the system at pH 8, the protonation states of Nδ1H for histidine (His), Oδ2(−) for aspartic acid (Asp), and Oε2(−) for glutamic acid (Glu) are adjusted accordingly. We first carry out simulations with Langevin dynamics using the CHARMM program with implicit solvent (EEF1), which is based on the CHARMM19 all-atom energy function with an effective Gaussian model to compute solvation free energy and nonbonding interactions.18,19 Second, the implicit-solvent REMD is performed using the Multiscale Modeling Tools for Structural Biology (MMTSB) toolset20 with a time step of 2 fs. The simulations are performed on 24 replicas within the 300−480 K temperature range with an exchange time step of 2 ps to allow the relaxation of the system. Each replica is simulated for a total of 40 ns, corresponding to a total simulation time of 0.96 μs. We further analyzed the structure obtained for the last 1000 exchanges of the production run in the trajectory at 300 K. The Kmeans clustering algorithm in the MMTSB toolset is used to cluster the structures based on mutual conformational similarity according to a root-mean-square deviation (RMSD) within 2 Å. Four representative structures with lowest energy and close to the center of clusters are identified, as shown in Figure 3. In order to study the ion effect, REMD with explicit solvent is applied to explore the ensemble of the Nvjp-1 protein structure in the presence of different concentrations of Zn2+ ions. These REMD simulations with explicit solvent and various ionic environments are carried out in NAMD 2.8.21 The CHARMM force field is applied22 to describe the interaction between atoms. The initial models used in the explicit-solvent REMD simulations are obtained from the representative structures obtained in the previous step with implicit solvent REMD. The entire atomistic structure is embedded completely in a TIP3 explicit water box with different ion concentrations (Zn2+/ protein = 0, 20, 50). A simulation time step of 2 fs is used, and the particle mesh Ewald function with a grid width of 1 Å is used to account for the electronic interactions. Conventional MD simulations were first performed on these systems in the NPT ensemble with a constant pressure of 1 atm at 300 K for 10 ns to equilibrate the structure in explicit solvent. The explicitsolvent REMD simulations were performed with the NVT ensemble at a temperature range of 300−480 K to simulate the ion effect on the protein structure. The computational cost is 768 CPUs for 96 replicas (8 CPUs for each replica) in XSEDE (Extreme Science and Engineering Discovery Environment). The equilibrium time within the exchange is 2 ps, with an average acceptance ratio of 36%. For each
ion concentration, every replica is simulated for a total of 2000 cycles until equilibrium is established. We then analyze the conformations from the trajectory at 300 K through the K-means clustering algorithm with 2 Å radius based on mutual similarity according to RMSD. The structures of selected clusters with the lowest energy are identified as the representative conformations. The secondary structure of each amino acid is identified using the STRIDE algorithm in Visual Molecular Dynamics (VMD) software.23 The SASA is calculated as the surface area that is accessible to a water molecule with a radius of 1.4 Å using VMD. Density Functional Theory Calculations. All histidine−metal complexes were fully optimized using DFT methods as implemented in the ORCA quantum chemistry package.24 The Becke threeparameter Lee−Yang−Parr functional (B3LYP)25,26 is adopted for optimizing the molecular structures, together with the 6-31G* basis set. To account for weak interactions, we include DFT-D3 dispersion correction with Becke−Johnson damping (D3BJ)27 and gCP-D3 energy correction,28 which has been recently developed. Additionally, the SMD model,29 a continuum solvation model based on the charge density of a solute molecule interacting with a continuum description of the solvent, is adopted. Elastic Network Model. The elastic network computational model is a bead−spring model with fcc lattice in which the equilibrated geometry and stiffness of the material are reflected by the spring constant (k) and the equilibrium length constant of the lattice structure (a). Both values are given as a function of the pH value in simulations. To be more specific, we use the fcc lattice to model the location of the mass beads (which locate at corners and face centers) and use the elastic spring to model the interaction between the nearest neighboring beads. The total energy function is given by E = ET =
∑ (r) φT
(1)
where ET is the deformation energy for stretching/compression of all the springs. Each energy term was defined by
φT(r ) = k(r − r0)2
(2)
where 2k is the spring stiffness, r0 = a/√2 is the equilibrium distance between two neighboring beads, and a is the lattice constant of the fcc lattice. By assuming the polymers have a compresion ratio of ν = 0.25, the Young’s modulus E of such a material composed of an fcc lattice and elastic springs is given by 1866
DOI: 10.1021/acsnano.6b07878 ACS Nano 2017, 11, 1858−1868
Article
ACS Nano k=
3aE 16
Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by grant number ACI-1053575.
(3)
where E is the Young’s modulus. According to our REMD simulation result that E = 2.5 GPa is obtained at pH = 5 and E = 10 GPa is obtained at pH = 8, we can simply use an approximate function as 2k = 2.88(pH − 5) + 6.5 (105 N/m) to implicitly model the effect of pH on spring stiffness and material stiffness. In the simulation of the coarsegrained model, we also implicitly model the effect of pH on material contraction by tuning the equilibrium length of the lattice constant as
a(pH) = −0.097 × pH + 1.18 (mm)
REFERENCES (1) Bryan, G. W.; Gibbs, P. E. Zinc - a Major Inorganic Component of Nereid Polychaete jaws. J. Mar. Biol. Assoc. U. K. 1979, 59, 969−973. (2) Lichtenegger, H. C.; Schoberl, T.; Ruokolainen, J. T.; Cross, J. O.; Heald, S. M.; Birkedal, H.; Waite, J. H.; Stucky, G. D. Zinc and Mechanical Prowess in the Jaws of Nereis, a Marine Worm. Proc. Natl. Acad. Sci. U. S. A. 2003, 100, 9144−9149. (3) Waite, J. H.; Broomell, C. C. Changing Environments and Structure−Property Relationships in Marine Biomaterials. J. Exp. Biol. 2012, 215, 873−883. (4) Broomell, C. C.; Mattoni, M. A.; Zok, F. W.; Waite, J. H. Critical Role of Zinc in Hardening of Nereis Jaws. J. Exp. Biol. 2006, 209, 3219−3225. (5) Birkedal, H.; Khan, R. K.; Slack, N.; Broomell, C.; Lichtenegger, H. C.; Zok, F.; Stucky, G. D.; Waite, J. H. Halogenated veneers: Protein Cross-Linking and Halogenation in the Jaws of Nereis, a Marine Polychaete Worm. ChemBioChem 2006, 7, 1392−1399. (6) Broomell, C. C.; Chase, S. F.; Laue, T.; Waite, J. H. Cutting Edge Structural Protein from the Jaws of Nereis Virens. Biomacromolecules 2008, 9, 1669−1677. (7) Broomell, C. C.; Zok, F. W.; Waite, J. H. Role of Transition Metals in Sclerotization of Biological Tissue. Acta Biomater. 2008, 4, 2045−2051. (8) Holten-Andersen, N.; Harrington, M. J.; Birkedal, H.; Lee, B. P.; Messersmith, P. B.; Lee, K. Y. C.; Waite, J. H. pH-induced MetalLigand Cross-Links Inspired by Mussel Yield Self-Healing Polymer Networks with Near-Covalent Elastic Moduli. Proc. Natl. Acad. Sci. U. S. A. 2011, 108, 2651−2655. (9) Fullenkamp, D. E.; He, L.; Barrett, D. G.; Burghardt, W. R.; Messersmith, P. B. Mussel-Inspired Histidine-Based Transient Network Metal Coordination Hydrogels. Macromolecules 2013, 46, 1167− 1174. (10) Srivastava, A.; Holten-Andersen, N.; Stucky, G. D.; Waite, J. H. Ragworm Jaw-Inspired Metal Ion Cross-Linking for Improved Mechanical Properties of Polymer Blends. Biomacromolecules 2008, 9, 2873−2880. (11) Brendel, V.; Bucher, P.; Nourbakhsh, I. R.; Blaisdell, B. E.; Karlin, S. Methods and Algorithms for Statistical Analysis of Protein Sequences. Proc. Natl. Acad. Sci. U. S. A. 1992, 89, 2002−2006. (12) Pappalardo, G.; Impellizzeri, G.; Bonomo, R. P.; Campagna, T.; Grasso, G.; Saita, M. G. Copper(ii) and Nickel(ii) Binding Modes in a Histidine-Containing Model Dodecapeptide. New J. Chem. 2002, 26, 593−600. (13) La Mendola, D.; Magrì, A.; Santoro, A. M.; Nicoletti, V. G.; Rizzarelli, E. Copper(II) Interaction with Peptide Fragments of Histidine−Proline-Rich Glycoprotein: Speciation, Stability and Binding Details. J. Inorg. Biochem. 2012, 111, 59−69. (14) Zhao, H.; Waite, J. H. Proteins in Load-Bearing Junctions: The Histidine-Rich Metal-Binding Protein of Mussel Byssus†,‡. Biochemistry 2006, 45, 14223−14231. (15) Gupta, R.; Dobritsa, S.; Stiles, C.; Essington, M.; Liu, Z.; Chen, C.-H.; Serpersu, E.; Mullin, B. Metallohistins: A New Class of Plant Metal-Binding Proteins. J. Protein Chem. 2002, 21, 529−536. (16) Hara, M.; Fujinaga, M.; Kuboi, T. Metal Binding by Citrus Dehydrin with Histidine-Rich Domains. J. Exp. Bot. 2005, 56, 2695− 2703. (17) Sugita, Y.; Okamoto, Y. Replica-Exchange Molecular Dynamics Method for Protein Folding. Chem. Phys. Lett. 1999, 314, 141−151. (18) Lazaridis, T.; Karplus, M. Effective energy function for proteins in solution. Proteins: Struct., Funct., Genet. 1999, 35, 133−152. (19) Brooks, B. R.; Brooks, C. L.; Mackerell, A. D.; Nilsson, L.; Petrella, R. J.; Roux, B.; Won, Y.; Archontis, G.; Bartels, C.; Boresch, S.; Caflisch, A.; Caves, L.; Cui, Q.; Dinner, A. R.; Feig, M.; Fischer, S.; Gao, J.; Hodoscek, M.; Im, W.; Kuczera, K.; et al. CHARMM: The Biomolecular Simulation Program. J. Comput. Chem. 2009, 30, 1545− 1614.
(4)
based on both REMD and experimental measurements, because the material length at pH = 8 shrinks to 40% of its length at pH = 2 and the initial lattice length is a = 0.5 mm for the model at pH = 7. Using the model, we can simulate the deformation of the composite with its surface exposed to different chemical conditions. It is noted that this model is a linearization around the bottom of a potential well in the energy/separation relationship, and it is important to note that this is strictly valid only for small perturbations around the equilibrium spacing (a). Expression and Purification of Nvjp-1. A codon-optimized Nvjp-1a1 isoform coding sequence was synthesized (Genescript) and cloned into a pET11 expression vector. The resulting pET11 Nvjp-1 construct was then transformed into BL21 DE3 cells. Expression cultures were inoculated and grown to an OD600 nm of 0.6 before being induced with 0.4 mM isopropyl β-D-1-thiogalactopyranoside in a shaking incubator for 12 h at 37 °C. Bacterial cells were pelleted, resuspended in 50 mM Tris-HCl, pH 8, 150 mM NaCl, and 100 μg/ mL lysozyme, and disrupted by sonication. After centrifugation, the insoluble pellet containing inclusion bodies was washed twice with 50 mM Tris-HCl, pH 8, 150 mM NaCl, and 1% Triton followed by one wash with 50 mM Tris-HCl, pH 8, and 150 mM NaCl. The pellet was resuspended in water, and the pH was adjusted to 3 with 1 N HCl followed by centrifugation. NaCl was then added to the supernatant to a final concentration of 100 mM, and the pH was adjusted to 8 with 1 N NaOH followed by centrifugation. The pellet was solvated in 6 mM HCl and dialyzed overnight against water. After lyophilization, Nvjp-1 was reconstituted to a final concentration of 60 mg/mL in water, centrifuged, and lyophilized a second time before being stored dry at −20 °C until use. Creation of Photo-Cross-Linked Nvjp-1 Hydrogels. Purified Nvjp-1 was reconstituted in water to a final concentration of 60 mg/ mL. A 15 μL amount of a 100 mM ammonium persulfate solution was added slowly to the reconstituted Nvjp-1 solution followed by the addition of tris(bipyridine)ruthenium(II) chloride (Ru(bpy)3Cl2) to a final concentration of 1 mM. This mixture was cast in a 2 × 0.5 × 0.15 cm3 polydimethylsiloxane mold and exposed to a strong white light source for 5 min. Ru(bpy)3Cl2 was removed by placing the photocross-linked hydrogels into a chelation buffer of 20 mM EDTA, 20 mM EGTA, and 1% β-mercaptoethanol, pH 8, for 24 h at room temperature. After removal of Ru(bpy)3Cl2, the hydrogels were swollen in 10 mM HCl prior to being treated with the respective salt solutions for 24 h at room temperature.
AUTHOR INFORMATION Corresponding Author
*E-mail:
[email protected]. ORCID
Rajesh R. Naik: 0000-0002-7677-928X Markus J. Buehler: 0000-0002-4173-9659 Notes
The authors declare no competing financial interest.
ACKNOWLEDGMENTS We acknowledge funding from the Air Force Office of Scientific Research (FA9550-11-1-0199). P.B.D. is adjunct faculty at Wright State University, Dayton, OH. This work used the NSF 1867
DOI: 10.1021/acsnano.6b07878 ACS Nano 2017, 11, 1858−1868
Article
ACS Nano (20) Feig, M.; Karanicolas, J.; Brooks, C. L. MMTSB Tool Set: Enhanced Sampling and Multiscale Modeling Methods for Applications in Structural Biology. J. Mol. Graphics Modell. 2004, 22, 377−395. (21) Phillips, J. C.; Braun, R.; Wang, W.; Gumbart, J.; Tajkhorshid, E.; Villa, E.; Chipot, C.; Skeel, R. D.; Kalé, L.; Schulten, K. Scalable Molecular Dynamics with NAMD. J. Comput. Chem. 2005, 26, 1781− 1802. (22) MacKerell, A. D.; Bashford, D.; Bellott; Dunbrack, R. L.; Evanseck, J. D.; Field, M. J.; Fischer, S.; Gao, J.; Guo, H.; Ha, S.; Joseph-McCarthy, D.; Kuchnir, L.; Kuczera, K.; Lau, F. T. K.; Mattos, C.; Michnick, S.; Ngo, T.; Nguyen, D. T.; Prodhom, B.; Reiher, W. E.; et al. All-Atom Empirical Potential for Molecular Modeling and Dynamics Studies of Proteins. J. Phys. Chem. B 1998, 102, 3586−3616. (23) Humphrey, W.; Dalke, A.; Schulten, K. VMD - Visual Molecular Dynamics. J. Mol. Graphics 1996, 14, 6. (24) Neese, F. The ORCA program system. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2012, 2, 73−78. (25) Becke, A. D. A New Mixing of Hartree−Fock and Local Density-Functional Theories. J. Chem. Phys. 1993, 98, 1372−1377. (26) Lee, C.; Yang, W.; Parr, R. G. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron Density. Phys. Rev. B: Condens. Matter Mater. Phys. 1988, 37, 785−789. (27) Grimme, S.; Ehrlich, S.; Goerigk, L. Effect of the Damping Function in Dispersion Corrected Density Functional Theory. J. Comput. Chem. 2011, 32, 1456−1465. (28) Kruse, H.; Grimme, S. A Geometrical Correction for the Interand Intra-Molecular Basis Set Superposition Error in Hartree-Fock and Density Functional Theory Calculations for Large Systems. J. Chem. Phys. 2012, 136, 154101. (29) Marenich, A. V.; Cramer, C. J.; Truhlar, D. G. Universal Solvation Model Based on Solute Electron Density and on a Continuum Model of the Solvent Defined by the Bulk Dielectric Constant and Atomic Surface Tensions. J. Phys. Chem. B 2009, 113, 6378−6396.
1868
DOI: 10.1021/acsnano.6b07878 ACS Nano 2017, 11, 1858−1868