Skip to main content
Advertisement
  • Loading metrics

Structural Refinement of Proteins by Restrained Molecular Dynamics Simulations with Non-interacting Molecular Fragments

  • Rong Shen,

    Affiliation Department of Biochemistry and Molecular Biology, Gordon Center for Integrative Science, University of Chicago, Chicago, Illinois, United States of America

  • Wei Han,

    Affiliation Beckman Institute, University of Illinois at Urbana–Champaign, Urbana, Illinois, United States of America

  • Giacomo Fiorin,

    Affiliation Institute for Computational Molecular Science, Temple University, Philadelphia, Pennsylvania, United States of America

  • Shahidul M. Islam,

    Affiliation Department of Biochemistry and Molecular Biology, Gordon Center for Integrative Science, University of Chicago, Chicago, Illinois, United States of America

  • Klaus Schulten,

    Affiliation Beckman Institute, University of Illinois at Urbana–Champaign, Urbana, Illinois, United States of America

  • Benoît Roux

    roux@uchicago.edu

    Affiliation Department of Biochemistry and Molecular Biology, Gordon Center for Integrative Science, University of Chicago, Chicago, Illinois, United States of America

Abstract

The knowledge of multiple conformational states is a prerequisite to understand the function of membrane transport proteins. Unfortunately, the determination of detailed atomic structures for all these functionally important conformational states with conventional high-resolution approaches is often difficult and unsuccessful. In some cases, biophysical and biochemical approaches can provide important complementary structural information that can be exploited with the help of advanced computational methods to derive structural models of specific conformational states. In particular, functional and spectroscopic measurements in combination with site-directed mutations constitute one important source of information to obtain these mixed-resolution structural models. A very common problem with this strategy, however, is the difficulty to simultaneously integrate all the information from multiple independent experiments involving different mutations or chemical labels to derive a unique structural model consistent with the data. To resolve this issue, a novel restrained molecular dynamics structural refinement method is developed to simultaneously incorporate multiple experimentally determined constraints (e.g., engineered metal bridges or spin-labels), each treated as an individual molecular fragment with all atomic details. The internal structure of each of the molecular fragments is treated realistically, while there is no interaction between different molecular fragments to avoid unphysical steric clashes. The information from all the molecular fragments is exploited simultaneously to constrain the backbone to refine a three-dimensional model of the conformational state of the protein. The method is illustrated by refining the structure of the voltage-sensing domain (VSD) of the Kv1.2 potassium channel in the resting state and by exploring the distance histograms between spin-labels attached to T4 lysozyme. The resulting VSD structures are in good agreement with the consensus model of the resting state VSD and the spin-spin distance histograms from ESR/DEER experiments on T4 lysozyme are accurately reproduced.

Author Summary

Knowledge of multiple conformational states of membrane transport proteins is a prerequisite to understand their function. However, the determination of atomic structures for all these states with conventional high-resolution approaches can be very challenging due to inherent difficulties in high yield purification of functional membrane transport proteins. Various complementary structural information of proteins in their native states can be obtained by a variety of biophysical and biochemical methods with site-directed mutations. Here, a novel restrained molecular dynamics structural refinement method is developed to help derive a structural model that is consistent with experimental data by incorporating all the experimental constraints simultaneously through the use of non-interacting all-atom molecular fragments. The method can be easily and effectively extended to incorporate many kinds of structural constraints from a variety of biophysical and biochemical experiments, and should be very useful in generating and refining models of proteins in specific functional states.

This is a PLOS Computational Biology Methods paper

Introduction

Membrane proteins can access multiple conformational states that are of functional importance. Knowledge of the three-dimensional (3D) structure of all the relevant conformations in atomic detail is necessary to understand how these membrane proteins accomplish their function at the molecular level. Despite great progress, it is often a challenge to determine high-resolution protein structures under physiological conditions using conventional approaches, e.g., x-ray crystallography [14], cryo-electron microscopy [5, 6] and nuclear magnetic resonance (NMR) spectroscopy [7, 8] due to the inherent difficulties in protein expression, purification and stability. The task is even more difficult when the goal is to determine the structure of specific functional states because such determination is possible only if one is able to stabilize the respective states through the control of experimental conditions or by introducing ligands or mutations. Alternatively, a wealth of complementary structural information of proteins in their native states can be obtained by biophysical and biochemical approaches. For example, electron paramagnetic resonance (EPR) spectroscopy can estimate solvent accessibilities of spin-labeled sites on a protein [912], double electron-electron resonance (DEER) spectroscopy, a powerful pulsed electron spin resonance (ESR) technique, can determine distance histograms between pairs of spin-labels covalently bound to a protein [13, 14], and solid-state NMR spectroscopy can investigate the orientation of transmembrane (TM) helices of membrane proteins and peptides in lipid environments [15, 16]. When combined with computational modeling methods and molecular dynamics (MD) simulations, the structural information from such experiments can be used for the determination and refinement of atomic structures of proteins in specific conformational states [1722]. However, translating experimental measurements into biologically meaningful and computationally efficient structural constraints is still quite challenging [20]. Residue-residue interactions derived from engineered salt or metal ion bridges are of particular interests as they can be translated directly into specific and simple spatial restraints, which can be incorporated into restrained MD simulations [23]. In particular, metal ion bridges (e.g., Cd2+, Zn2+ and Mg2+) can be easily formed between pairs of cysteines, histidines, as well as glutamates and aspartates when they are within atomic proximity [24, 25].

Studies of voltage-sensing domains (VSDs) in various ion channels provide a good illustration of the benefits from using this type of strategy [2632], as salt or metal ion bridges can break and reform during the functional cycle of a working channel, avoiding trapping the protein in non-functional or distorted conformations due to their low strength [32]. The VSD structure is formed by a small bundle of four antiparallel TM helices, S1 to S4. A series of basic residues distributed along the helix S4 provide the “gating charge” that senses the changes in membrane potential. Upon depolarization, the VSDs rearrange from resting state into active state, triggering the opening of the central pore domain of the channel [33, 34]. The scenario that conveys most accurately the voltage-activated conformational change of the VSD is the “sliding helix” mechanism in which the S4 segment translates principally along its main axis while maintaining its helical conformations [3537]. It is believed that during the voltage-activated transition, the basic residues along S4 form sequential salt bridge interactions with acidic residues located along S1–S3.

The x-ray crystal structures of the VSD in the active state have been determined for a few voltage-gated K+ [38, 39] and Na+ [4042] channels. However, there is currently no atomic structure of the VSD of voltage-gated ion channels in the resting state. The only available crystal structure of the resting state VSD is that of the voltage-sensitive phosphatase from C. intestinalis (Ci-VSP), which is a distant homolog of the VSD of voltage-gated ion channels [4]. Using structural restraints derived from cysteine-cysteine Cd2+ bridges, Henrion et al. constructed one open and four different closed structures of the VSD of the Shaker potassium channel by Rosetta modeling and MD simulations [32]. Vargas et al. built and simulated several models of the VSD of the Kv1.2 channel in the resting state with experimentally determined salt and metal ion bridges, and the results of these restrained MD simulations were then used to generate a consensus model of the resting state conformation of the VSD, which is consistent with a wide range of available experimental data [23, 43] and displays a high degree of similarity to the atomic models of the resting state VSD constructed by different research groups with different methods [4446]. In these two studies, the experimentally derived structural constraints, e.g., metal ion bridges, were imposed explicitly by performing proper site-directed mutations and adding specific bridging ions, instead of being represented as simplified interatomic distance restraints [34, 46]. However, each of these structural constraints was applied individually in independent models of the VSD. It was not possible to explicitly model all the structural constraints at once as one site of the protein might be involved in different engineered structural constraints and would cause unphysical steric clashes.

In the present study, we develop a novel structural refinement method by introducing the concept of non-interacting molecular fragments, which are used in the context of restrained MD simulations that can integrate the information from multiple experimental constraints simultaneously. The molecular fragments (e.g., metal ion bridges or spin-labels) are defined with all atomic details (Fig 1A). The fragments are anchored to the protein via their overlapping backbone atoms (N, C, O and Cα) using harmonic restraints (Fig 1B). The interaction within each molecular fragment is treated realistically, while there is no interaction between different molecular fragments. To facilitate the refinement of atomic models of proteins, the molecular fragment method has been implemented in NAMD [47] that is efficient in MD simulations, and has been assisted further by a VMD [48] plugin designed to provide a user-friendly interface for constructing the models (S1 Fig). To test and illustrate the feasibility of the method, we refined the models of the resting state VSD of the Kv1.2 channel from Khalili-Araghi et al. [44] in an explicit lipid bilayer with multiple molecular fragments based on engineered salt and metal ion bridges between the four TM helices (S1–S4) of the VSD [2832]. The VSD provides a good illustration for the present method because the 3D structure of the active state is known from crystallography, and experimental data must be used to determine the conformation of the resting state (Fig 1C). The results of the restrained MD simulations show that the structural restraints can be satisfied simultaneously without large conformational changes of the VSD from its starting configuration. The final structures are in good agreement with each other and with a previously established consensus model [23, 43], and are stable in subsequent simulations without the restraints. The non-interacting molecular fragment method has also been extended to incorporate the mean-field restrained-ensemble MD (reMD) simulation method with methanethiosulfonate spin-labels (MTSSL) as molecular fragments and the ESR/DEER distance histograms as structural constraints [4951]. A total of 51 pairs of experimental distance histogram restraints were applied to 34 explicit spin-labels inserted at different positions in a single T4 lysozyme (T4L), and each spin-label was expanded into 25 copies. It is found that all the experimental distance histograms can be accurately reproduced simultaneously in the reMD simulations.

thumbnail
Fig 1. Schematic of the molecular fragments method.

(A) Each structural constraint (e.g., a metal ion bridge or a spin-label) is present as a molecular fragment with all atomic details in the system. (B) The residue(s) in the molecular fragment (colored in red) is/are attached to the targeting residue(s) in the wild-type protein, via harmonic restraints, with their backbone atoms (N, C, O and Cα) staying on top of each other, respectively. The residue in the molecular fragment does not have interactions with its targeting residue in the wild-type protein (residue i), and the interactions between the backbone atoms of the residue in the molecular fragment and the backbone atoms of the two nearby residues in the wild-type protein (residues i-1 and i+1; embraced by red dashed lines) are also turned off. (C) During the restrained MD simulations, the interactions within each molecular fragment are evaluated accurately triggering the conformational changes of the wild-type protein, while different molecular fragments do not have interactions with one another.

https://doi.org/10.1371/journal.pcbi.1004368.g001

Results and Discussion

Four models of the resting state VSD of Kv1.2 with different structural constraints were generated (Table 1). In the first model (model-1), shown in Fig 2A, we simultaneously applied the first five structural constraints between S1 and S4, S2 and S4, and S2 and S3. Vargas et al. separately imposed these five structural constraints in their refinement of the resting state VSD of the Kv1.2 channel [23]. It was shown that each of the restraints could be accommodated by the VSD with little rearrangement of the backbone of S1–S4, and the requirement of the five interactions could all be satisfied in their final consensus structural model. As a result, we included the first five structural constraints in all of the four models. The last two structural constraints were reported recently [32], thus they were further incorporated separately or all at once in model-2 to model-4, in combination with the first five structural constraints, as a comparison.

thumbnail
Table 1. Details of the models of the VSD of the Kv1.2 channel.

https://doi.org/10.1371/journal.pcbi.1004368.t001

thumbnail
Fig 2. Starting configurations of the four models of the VSD with different molecular fragments.

(A) model-1. (B) model-2. (C) model-3. (D) model-4. The four TM helices, S1 (gray), S2 (yellow), S3 (red) and S4 (blue), are represented in transparent ribbons. Molecular fragments, I177C–Cd2+–R294C, I230C–Cd2+–R294C, I230D–Mg2+–F267D, F233W/E236–R294K, I230H–Zn2+–A291H, I268C–Cd2+–A287C, and T269C–Cd2+–A291C, are colored in red, yellow, green, magenta, cyan, orange, and blue, respectively, with amino acids being represented in sticks and ions in spheres.

https://doi.org/10.1371/journal.pcbi.1004368.g002

Practically, the heavy backbone atoms of the residue in the molecular fragment (colored in red in Fig 1B) are anchored to the respective heavy backbone atoms of the assigned residue in the VSD (residue i in Fig 1B). The interaction between the residue of the molecular fragment and its targeting residue (residue i) of the VSD, and the interaction between the backbone atoms of the molecular fragment residue and the backbone atoms of the two nearby residues of residue i (residue i-1 and residue i+1) are turned off. While the interactions within each molecular fragment are unaltered and treated realistically, the molecular fragments of the different constraints do not see one another (i.e., they are “non-interacting”) to avoid unphysical steric clashes.

In the initial structure of model-1, the distance between the Cβ atoms of the two cysteine residues in the molecular fragment MF-1 (I177C–Cd2+–R294C) is 9.65 Å and the distances between the Cd2+ ion and the sulfur atoms of the cysteine residues are 5.10 Å. The corresponding Cβ-Cβ and Sγ-Cd2+ distances are 10.26 and 3.69 Å, respectively, in the molecular fragment MF-2 (I230C–Cd2+–R294C). The distances are too large for a Cys–Cd2+–Cys bridge, as the geometric restraint of a Cys–Cd2+–Cys bridge requires that the Cβ atoms of the two cysteines be apart from each other by not more than 9 Å [32].

To improve the model, harmonic distance restraints between the Cd2+ ion and the Sγ atoms, centered at 2.6 Å [24], and harmonic angle restraints between Sγ-Cd2+-Sγ, centered at 180°, were applied in the restrained MD simulations containing cysteine-cysteine Cd2+ bridges (Table 2). To allow the VSD and the molecular fragments to adjust their conformations smoothly, the force constants for the harmonic restraints were gradually turned on during the simulations (S1 Table). The average distances between Sγ and Cd2+ and between Cβ and Cβ for MF-1 had been reduced to 2.40 ± 0.05 Å and 6.58 ± 0.19 Å, respectively, prior to the last 10 ns of the simulations. For MF-2, the average Sγ-Cd2+ distance was 2.52 ± 0.07 Å, which is similar to that for MF-1, while the average Cβ-Cβ distance was 8.24 ± 0.14 Å, 1.66 Å longer than that for MF-1. It is expected that the Cβ-Cβ distance exhibits a wider spread than the ligand-metal (e.g., Sγ-Cd2+, Oδ-Mg2+ and Nε-Zn2+) and salt bridge (e.g., Nζ-Oε and Nζ-Oδ) distances due to the multiple rotameric states of the side chains of the interacting residues (Figs 3 and 4).

thumbnail
Table 2. Force field parameters for the metal ion bridges.

https://doi.org/10.1371/journal.pcbi.1004368.t002

thumbnail
Fig 3. Final configurations of the four models of the VSD with different molecular fragments.

The final configuration (left panel) and the correlation between the Cβ-Cβ distance and the Sγ-Cd2+, Oδ2-Mg2+, Nε-Zn2+, or Nζ-Oε2 distance in each molecular fragment (right panel) of (A) model-1, (B) model-2, (C) model-3, and (D) model-4. The last 10 ns trajectories of the restrained MD simulations were used to calculate the average values and standard deviations of the distances.

https://doi.org/10.1371/journal.pcbi.1004368.g003

thumbnail
Fig 4. The final coordinates of the molecular fragments in model-3 after the restrained MD simulation.

(A–F) The configurations of the four TM helices before and after the restrained MD simulation are shown in transparent and solid ribbons, respectively.

https://doi.org/10.1371/journal.pcbi.1004368.g004

In addition to the Cd2+ ion bridge, two other metal (Mg2+ and Zn2+) ion bridges were also considered in the models. In the molecular fragment MF-3 (I230D–Mg2+–F267D), the Cβ atoms of the aspartate residues are initially 12.59 Å apart. The distances between the Mg2+ ion and the Oδ2 atoms and the Oδ2-Mg2+-Oδ2 angle were harmonically restrained during the simulations, with an equilibrium bond distance of 2.1 Å [25] and an equilibrium angle of 180°, to form the aspartate-aspartate Mg2+ bridge (Table 2). The average Cβ-Cβ and Mg2+-Oδ2 distances as seen in the last 10 ns of the trajectory were 7.51 ± 0.08 Å and 1.94 ± 0.07 Å, respectively. It is also seen that the Mg2+-Oδ2 distance is very stable in all of the four models (Fig 3), attributed to the strong bidentating interactions between Mg2+ and carboxylate groups of the aspartate residues [25]. For the histidine-histidine-Zn2+ bridge in the molecular fragment MF-5 (I230H–Zn2+–A291H), the imidazole nitrogen atom Nε was assumed to coordinate the Zn2+ ion, as it appears to be the configuration with the highest propensity in the protein data base [24]. Furthermore, so-called improper torsion and dihedral restraints were applied to restrain the Zn2+ ion within the planes of the imidazole rings (Table 2). Harmonic restraints were applied to control the Cε1-Nε-Zn2+ and Cδ2-Nε-Zn2+ angles, employing geometrical characteristics according to the crystal structures (e.g., PDB IDs: 1CBX, 1I4P, 1STE and 3CPA) containing the histidine-histidine-Zn2+ bridge [5255]. Noting that the Zn2+ ion could also be bound by the Nδ atom, as observed in the crystal structures, we tried to restrain the Zn2+-Nδ bonds in our models and found that the Zn2+-to-Nδ distance could easily be satisfied, but that the orientation of the Zn2+-Nδ bonds with respect to the imidazole rings was hard to be maintained.

It has been shown that the double mutant F290W/R362K can further stabilize the resting state of the Shaker channel than the single mutant F290W [30]; apparently, the double mutant enables the formation of an electrostatic interaction between R362K and E293, which is one helix turn below F290W, in the resting state. To incorporate this lysine-glutamate salt bridge into the models of the VSD resting state of the Kv1.2 channel, we constructed a molecule fragment of F233W–R294K (MF-4). The amino group of the R294K side chain and the carboxylate group of the E236 side chain are initially far away from each other (further than 16 Å). To form the salt bridge, a harmonic restraint was applied between the Nζ atom of R294K and the Oε2 atom of E236 with an equilibrium bond distance of 3.5 Å. The Nζ atom moved towards the Oε2 atom quickly during the simulation and the average Nζ-Oε2 distance reduced to 3.72 ± 018 Å for the final 10 ns, with a decrease of the distance between the Cβ atoms of E236 and R294K within roughly 5 Å. Superposition of the initial and final configurations of the VSD of model-1 using the STAMP structural alignment algorithm [56] reveals little rearrangement of the protein (Fig 5A). The backbone root-mean-square deviations (RMSDs) between the refined structure of model-1 with its initial structure and the consensus model, in the TM regions of the S1–S4 helices, are 2.2 Å and 3.1 Å (Table 3), respectively, indicating that the five residue-residue interactions applied by incorporating the proper molecular fragments into the model can be fulfilled simultaneously without large conformational changes of the protein backbone.

thumbnail
Fig 5. Superimposition of the initial and final configurations of the four models of the VSD.

(A) model-1, (B) model-2, (C) model-3, and (D) model-4. The initial and final configurations of the backbone atoms are represented in transparent and solid sticks, respectively.

https://doi.org/10.1371/journal.pcbi.1004368.g005

thumbnail
Table 3. The backbone RMSD values of the TM helices S1–S4 between different models.

https://doi.org/10.1371/journal.pcbi.1004368.t003

Two additional interactions between residues in the S3 and S4 helices, which have been recently found to arise in the resting state VSD of the Shaker channel [32], were further considered in our models. We incorporated another Cd2+ ion bridge molecular fragment MF-6 (I268C–Cd2+–A287C) and MF-7 (T269C–Cd2+–A291C) into model-2 and mode-3 (Fig 2), respectively, together with the initial five molecular fragments. As shown in Fig 3B, the average Sγ-Cd2+ and Cβ-Cβ distances of MF-6 had been reduced to 2.46 ± 0.06 Å and 6.34 ± 0.18 Å from 8.52 Å and 15.02 Å, respectively, for the final 10 ns. However, the average Sγ-Cd2+ distances of MF-1 and MF-2 in model-2 are slightly larger (by ~ 0.5 Å) than the distance in model-1, with increased standard deviations, and the average Cβ-Cβ distance of MF-4 is about 1.5 Å larger in model-2 than in model-1. The distance variations are consistent with the conformational rearrangement of the model. Comparison of the starting and final conformation of model-2 reveals an inward movement of S4, and an outward translation and a clockwise rotation (top view) of S3, while the S1 and S2 helices mainly retain their conformation (Fig 5B). The inward movement of S4, which resulted from the newly applied structural restraint MF-6, changed the configuration of the molecular fragments MF-1, MF-2 and MF-4, restraining the relative positions between S1 and S4, and S2 and S4. Similar conformational changes were observed for model-3 (Fig 5C), using the molecular fragment MF-7 to regulate the interaction between residues in S3 and S4. However, the S3 helix moved further outward and the S4 helix moved less inward in model-3 than in model-2, and the residue-residue interactions in the molecular fragments were well fulfilled (Figs 3C and 4). The two structural restraints (MF-6 and MF-7) were both applied in our last model, model-4, (Fig 2D). It was found that the S3 helix kinked in the center at the end of the simulation, although the structural restraints of MF-6 and MF-7 were satisfied (Figs 3D and 5D). We propose that the resides involved in the last two molecular fragments connecting the S3 and S4 helices are too close to each other, and the S3 helix is not able to adjust its conformation to satisfy the two structural constraints at the same time while maintaining an alpha-helical structure. However, they can be satisfied separately with slight conformational changes of the VSD as shown in model-2 and model-3. The RMSDs of the backbone atoms in the TM region of the S1–S4 helices between the four refined models and the models of Khalili-Araghi et al. [44], Vargas et al. [23] and Jensen et al. [45] lie within 1.7–3.8 Å of one another (Table 3), revealing a high degree of similarity between these structures at the atomic level [43]. In addition, the conformation of the four refined models remain stable during an unbiased MD simulation of 100 ns without the structural restraints (see S2 Fig), implying that the structural constraints will not displace the conformation of the VSD away from its “native” state(s). The Protein Data Bank (PDB) coordinates of the TM region of the refined model-3 are provided in Supplementary Material.

The present methodology with multiple non-interacting molecular fragments can also be used in the treatment of DEER data. DEER spectroscopy is a powerful electron spin resonance technique allowing the determination of distance histograms between pairs of spin-labels attached to a protein in their native environment. Recently, a mean-field restrained-ensemble molecular dynamics (reMD) simulation method was introduced to exploit ESR/DEER distance histograms in protein structural refinement [4951]. The method was built on the framework of the multiple-copy locally enhanced sampling method, in which each spin-label side chain was expanded into N copies, yielding a total of N2 spin-spin distances for each pair of spin-labels. Energy restraints were then implemented to couple the ensemble of copies of the spin-labels and enforce the experimental ESR/DEER distance histograms. In the reMD simulations, the different copies of a given spin-label do not have interactions with one another while the copies from different spin-labels can interact with each other. It requires that the spin-labels should be kept apart with certain distances from one another, downgrading the efficiency of the reMD simulation method. For example, up to 51 pairs of spin-spin distance histograms have been measured for a small benchmark protein of 165 residues, T4L, involving a total of 34 different site-directed spin-labeling sites. To apply all of the 51 experimental restraints simultaneously, four separate copies of the protein T4L were included in the simulation system to prevent steric clashes between the spin-labels [50].

The non-interacting molecular fragment method makes it possible to include all of the spin-labels together in a single copy of the protein in reMD simulations. As shown in Fig 6A, spin-labels were introduced at 34 positions in T4L, and the spin-label at each site has 25 copies. The 51 distance distribution data from DEER were then simultaneously used in the reMD simulation [49], which enforce an energy restraint to match the calculated distance distribution to agree with those from DEER (Fig 6B). All the 51 distance distributions data are shown in the Supplementary Material, S3 and S4 Figs. It should be noted that the unbiased reMD simulation without any restraint was not able to converge accurately toward the experimental distance histograms, especially for those histograms having complex structures (S3 and S4 Figs), as discussed in the previous studies [49, 50]. It is worth mentioning that the present molecular fragment method can also be used in a simpler reMD treatment of DEER data in which the spin-labels are represented by dummy spin-labels [4951]. The reMD simulation of the dummy spin-label attached to T4L was also found to rapidly match the distance distribution with those obtained from EPR/DEER within a few ns (data not shown).

thumbnail
Fig 6. Configurations of the spin-labels attached to T4 lysozyme.

(A) A snapshot of the T4 lysozyme system with 34 MTSSL spin-labels. The spin-label at each site has 25 copies and highlighted in the ball-and-stick representation. (B) Distance histograms of four pairs of spin-labels from ESR/DEER experiments (black solid line), and the molecular dynamics simulations without (blue dashed line) and with (red dashed line) an energy restraint. Please see the Supplementary Material, S3 and S4 Figs for the other 47 distance distributions.

https://doi.org/10.1371/journal.pcbi.1004368.g006

In summary, the atomic resolution structures of membrane transport proteins in different states are critical for developing an understanding of the proteins’ biological function at the molecular level. However, the determination of the high-resolution structures using conventional structure determination methods encounters different types of challenges and is often not feasible, especially for membrane proteins. We have developed a novel method based on non-interacting molecular fragments with restrained molecular dynamics to refine protein structures using information from residue-residue interactions derived from a wide range of biophysical and biochemical experiments. The multiple independent structural constraints are incorporated simultaneously on the basis of individual non-interacting molecular fragments represented in atomic detail. Using the new method, we have refined the atomic resolution structure of the VSD of the Kv1.2 channel in the resting state based on different sets of experimental restraints. A number of models were generated using different subset of experimental constraints. The method could efficiently generate a 3D structural model of the VSD that satisfies all the imposed experimental constraints simultaneously. The final refined models are quite similar to the Vargas et al. consensus model [23, 43] with the backbone RMSD between the final models and the consensus model, in the TM region of the S1–S4 helices, seen to lie between 2.6 and 3.1 Å. Extensions and applications of the molecular fragment method to incorporate other advanced simulation methods and structural constraints can also be achieved easily and conveniently. An illustrative test of the spin-spin distance distribution of T4L shows that the molecular fragment method is compatible with the mean-field reMD simulation method by treating each copy of the spin-labels as an individual molecular fragment. The resulting distance histograms of 51 pairs of spin-labels from a single copy of the protein accurately matched the experimental data simultaneously by applying a novel energy restraint to the ensemble of spin-labels.

Materials and Methods

Molecular fragment method

The molecular fragment method refines the VSD structural models by applying in the same simulations different structural constraints derived from individual experiments [2832]. Each structural constraint corresponds to a molecular fragment containing a salt bridge or a metal ion bridge that connects two remote residues of the VSD. These molecular fragments were represented with all atomic details (Fig 1) and attached to assigned residues of the VSD by keeping the backbone atoms (N, C, O and Cα) of the residues in the molecular fragments and their targeting residues in the VSD staying on top of each other through harmonic restraints centered at 0 Å. The bridging metal ions were placed at the centers of the fragments, coordinating to the Sγ atom of cysteine, the Nε atom of histidine, the Oδ2 atom of aspartate or the Oε2 atom of glutamate depending on types of fragments (S1 Movie).

For the sake of simplicity, the cysteine, histidine, aspartate and glutamate residues involved metal ion bridges were assigned in the deprotonated state carrying a negative charge [23], such that each molecular fragment is electrically neutral. We have built a VMD plugin to introduce the molecular fragments to the model and generate input files for the subsequent MD simulations (S1 Fig). This plugin is available in the Supplementary Material and at (www.ks.uiuc.edu/Research/MolFragReMD/).

All simulations containing molecular fragments were performed with NAMD [47]. To avoid steric clashes when simulating multiple molecular fragments, NAMD was modified such that any interaction between the fragments themselves was removed (“switched off”) during simulations. As heavy backbone atoms of the fragments were forced to stay on the top of the respective atoms of the targeting residues of the protein, the interactions between the overlapping parts of the fragments and the protein atoms were also removed (Fig 1). To further avoid any potential bad contact, an effective distance, , was introduced to NAMD to calculate the interaction of pairs of atoms between the molecular fragments and the protein, where Δ is an input parameter chosen by the user. With this effective distance, unphysical small distances between pairs of atoms are allowed to occur without introducing abrupt fluctuations in the atomic forces that destabilize the simulations.

A technical issue associated with molecular fragment simulations is how to eliminate electrostatic interactions between different fragments. When using the particle mesh Ewald (PME) method [57], long-range electrostatic forces are calculated in reciprocal space. If the partial charges of all the fragments were all accumulated on the PME grid, the long-range electrostatic interactions between two fragments would be difficult to subtract in a pairwise manner. Because the fragments are restrained around the target residues of the protein, the latter are instead a good approximation for calculating long-range interactions. In this approach, only the protein’s charges are accumulated on the PME grid when calculating long-range electrostatic forces, while the fragments interact with their surroundings with truncated Coulomb and van der Waals potentials.

Construction of the starting models

The structure of a single VSD (residues: 158–308) of the Kv1.2 channel in the resting state proposed by Khalili-Araghi et al. [44], subsequently assembled from the atomic model of Pathak et al. [22], was used to generate our starting models with VMD [48]. Here, seven non-interacting molecular fragments (labeled MF-1 to MF-7) were considered to model the resting state of the VSD: MF-1 is a cysteine-cysteine Cd2+ bridge involving residues I241C (S1) and R362C (S4) in Shaker [28]; corresponding residues in the Kv1.2 channel are I177C (S1) and R294C (S4). MF-2 is a cysteine-cysteine Cd2+ bridge involving residues I287C (S2) and R362C (S4) in Shaker [28]; corresponding residue in Kv1.2 are I230C (S2) and R294C (S4). MF-3 is an aspartate-aspartate Mg2+ bridge involving residues I287D (S2) and F324D (S3) [29], equivalent to I230D (S2) and F267D (S3) in the Kv1.2 channel. MF-4 is a lysine-glutamate salt bridge involving F290W (S2), R362K (S4) and E293 (S2) in Shaker [30]; corresponding residue in Kv1.2 are F233W (S2), R294K (S4) and E236 (S2). MF-5 is a histidine-histidine Zn2+ bridge involving residues I287H (S2) and A359H (S4) in Shaker [31]; corresponding residues in Kv1.2 are I230H (S2) and A291H (S4). MF-6 is a cysteine-cysteine Cd2+ bridge involving residues I325C (S3) and A355C (S4) in Shaker [32]; corresponding residue in Kv1.2 are I268C (S3) and A287C (S4) in the Kv1.2 channel. MF-7 is a cysteine-cysteine Cd2+ bridge involving residues T326C (S3) and A359C (S4) in Shaker [32]; corresponding residue in Kv1.2 are T269C (S3) and A291C (S4). In the case of MF-4, the double mutant F290W/R362K supposedly enables an electrostatic interaction between R362K and E293 in the Shaker channel. Four models of the VSD in the resting state were generated by considering various combinations of these interactions (Table 1). The histidine residues were modeled in a protonation state with a positive charge, and all other ionizable residues were assigned in their default ionization state, i.e., Glu and Asp carrying a negative charge, Lys and Arg carrying a positive charge, as in the atomic models of Khalili-Araghi et al [44]. The VSD was then embedded into a pre-hydrated 1-palmitoyl-2-oleoyl-sn-glycero- 3-phosphocholine (POPC) lipid bilayer. The membrane normal axis was aligned along the z-axis with the center of the membrane at z = 0. An additional sodium ion was inserted into the bulk water to make the system in an electrically neutral state.

To simulate DEER data, the x-ray crystal structure of T4L (PDB ID: 2LZM) was used to generate the starting T4L system with VMD [48]. The protein was solvated into a cubic water box with 8 chloride ions. The ionizable residues of T4L were assigned in their default ionization state and the final system is electrically neutral.

Restrained molecular dynamics simulations

The protein, phospholipids and ions were described with the CHARMM36 force field [5860], and water molecules with the TIP3P model [61]. Harmonic restraints with soft force constants were applied to maintain the conformation of the metal bridges [24, 25]. Harmonic restraints were also applied to the heavy backbone atoms (N, Cα, C and O) of the molecular fragments and the corresponding atoms of their targeting residues (Table 2).

Following 5000 steps of energy minimization, each system was initially equilibrated for 1 ns with the protein and the molecular fragments being restrained to fill the gap between the protein and the lipids. The systems were then simulated for 20 ns with the harmonic restraints applied to the metal bridges being gradually increased. Secondary structure restraints were also applied to the four TM helices of the VSD to avoid large distortions of the protein at the initial stages of the refinement, and the force constant was gradually decreased. Another 80 ns simulation was carried out for each system with constant restraints to test the stability of the refined VSD and the molecular fragments (S1 Table). The refined systems were further simulated for 100 ns without the structural constraints, and the backbone atoms of the protein were restrained during the first 10 ns using a force constant of 1 kcal/mol/Å2.

A time step of 2 fs was used with the bonds involving hydrogen atoms being constrained using the SHAKE algorithm [62]. Electrostatic interactions were calculated using the PME method [57] with a grid density of at least 1 bin per Å3, and the van der Waals interactions were smoothly switched off after 10 Å and truncated beyond 12 Å. Periodic boundary conditions were imposed in all directions. The temperature of the systems was controlled at 300 K using the Langevin dynamics and the pressure was kept at 1 atm using the Nose-Hoover Langevin piston method [63, 64].

Restrained-ensemble molecular dynamics (reMD) simulations

Spin-labels, consisting of disulfide bond linked MTSSL (1-oxyl-2,2,5,5-tetramethylpyrroline-3-methyl methanethiosulfonate) and cysteine, were attached at 34 sites of T4L using the VMD plugin we built. The spin-label at each site has 25 copies, and each of the copies was assigned as a separate molecule fragment. The heavy backbone atoms (N, C, O and Cα) of the spin-labels and those of the targeting residues of T4L were forced to stay on top of each other, respectively, through harmonic restraints with a force constant of 10 kcal/mol/Å2 and an equilibrium distance of 0 Å. The interactions between the spin-labels and the rest of the system were scaled by 1/25. Because of the coding design of NAMD, one can only define at most 255 non-interacting molecular fragments, while there are a total of 850 copies (34×25) of spin-labels in the T4L system. To circumvent this constraint, several spin-labels that were separated by a distance larger than the non-bonded cutoff were grouped in the same molecular fragment. For example, the s-th copy of spin-label at site 59 and the s’-th copy of spin-label at site 128 could be classified into the same molecular fragment, while the s-th and s’-th copies of spin-label at site 59 or 128 were classified into two different molecular fragments, thus each copy of the spin-labels has no interactions with one another.

The distance histogram for each pair of spin-labels was restrained toward the experimental data from ESR/DEER using an ensemble energy restraint introduced by Roux and coworkers [49, 50] where and Hij(n) are the n-th bin of the histogram for the ij pair of spin-labels from the reMD simulations and experiments, respectively. A smooth differentiable Gaussian of width σ was used to construct the histograms where N is the number of copies of a given spin-label, Δr is the bin width, is the distance between the oxygen atoms of the s-th copy of spin-lable i and the s’-th copy of spin-lable j. The histograms from the simulations and experiments were properly normalized, namely, for all the ij pair of spin-labels. The energy restraint has been introduced into NAMD under the collective variables (colvars) module. The force constant K was set to 500 kcal/mol/Å2, and the natural spread σ of the Gaussian was set to 1.1 Å. The bin width Δr was set to 1 Å and n = 1, …, 60. The backbone of the T4L was harmonically restrained with a force constant of 100 kcal/mol/Å2 to prevent any large displacement and distortion.

The system was energy minimized and equilibrated for 2 ns without the energy restraint, followed by a 1 ns production run with the energy restraint. The simulations were performed in the NPT ensemble with a temperature of 300 K and a pressure of 1 atm using the Nose-Hoover Langevin piston method [63, 64]. The time step of the simulations was set to 1 fs and the bonds involving hydrogen atoms were constrained with the SHAKE algorithm [62]. The PME method was used to calculate the electrostatic interactions [57], and the van der Waals interactions were smoothly switched off from 10 to 12 Å. For the sake of completeness, the reMD simulation was also carried out using dummy spin-labels according to the previously established protocol [4951].

All the features stated above, including the molecular fragment method and the distance histogram restraints, have been implemented in a modified version of NAMD. The source code and a brief tutorial of the molecular fragment method are available at (www.ks.uiuc.edu/Research/MolFragReMD/).

Supporting Information

S1 Fig. The interface of the VMD plugin used for constructing models with multiple molecular fragments.

https://doi.org/10.1371/journal.pcbi.1004368.s001

(TIF)

S2 Fig. The trajectories of the backbone RMSDs of the TM region of the four VSD models as a function of simulation time and the superimposition of the starting and final configurations of the models in the extended MD simulations without the structural constraints.

https://doi.org/10.1371/journal.pcbi.1004368.s002

(TIF)

S3 Fig. The distance histograms (h1–h28) for the spin-labels attached to T4 lysozyme from the ESR/DEER experiments and the mean-field restrained-ensemble simulations.

https://doi.org/10.1371/journal.pcbi.1004368.s003

(TIF)

S4 Fig. The distance histograms (h29–h51) for the spin-labels attached to T4 lysozyme from the ESR/DEER experiments and the mean-field restrained-ensemble simulations.

https://doi.org/10.1371/journal.pcbi.1004368.s004

(TIF)

S1 Table. Restraints applied to the metal bridges and the four transmembrane helices.

https://doi.org/10.1371/journal.pcbi.1004368.s005

(PDF)

S1 Movie. Animation showing the process of constructing models with molecular fragments.

https://doi.org/10.1371/journal.pcbi.1004368.s006

(MPG)

S1 PDB. Backbone of the S1–S4 helices of the refined model-3 for the voltage-sensing domain of Kv1.2.

https://doi.org/10.1371/journal.pcbi.1004368.s007

(PDB)

S1 VMD-Plugin. Scripts and usage of the VMD plugin used for constructing models with multiple molecular fragments.

https://doi.org/10.1371/journal.pcbi.1004368.s008

(ZIP)

Acknowledgments

RS would like to thank Ernesto Vargas and Fatemeh Khalili-Araghi for helpful discussions.

Author Contributions

Conceived and designed the experiments: RS WH GF SMI KS BR. Performed the experiments: RS WH GF SMI. Analyzed the data: RS WH GF SMI. Contributed reagents/materials/analysis tools: RS WH GF SMI. Wrote the paper: RS WH GF SMI KS BR.

References

  1. 1. Cuello LG, Jogini V, Cortes DM, Pan AC, Gagnon DG, Dalmas O, et al. Structural basis for the coupling between activation and inactivation gates in K+ channels. Nature. 2010;466: 272–275. pmid:20613845
  2. 2. Doyle DA, Morais Cabral J, Pfuetzner RA, Kuo A, Gulbis JM, Cohen SL, et al. The structure of the potassium channel: molecular basis of K+ conduction and selectivity. Science. 1998;280: 69–77. pmid:9525859
  3. 3. Cuello LG, Jogini V, Cortes DM, Perozo E. Structural mechanism of C-type inactivation in K+ channels. Nature. 2010;466: 203–208. pmid:20613835
  4. 4. Li Q, Wanderling S, Paduch M, Medovoy D, Singharoy A, McGreevy R, et al. Structural mechanism of voltage-dependent gating in an isolated voltage-sensing domain. Nat Struct Mol Biol. 2014;21: 244–252. pmid:24487958
  5. 5. Cao E, Liao M, Cheng Y, Julius D. TRPV1 structures in distinct conformations reveal activation mechanisms. Nature. 2013;504: 113–118. pmid:24305161
  6. 6. Liao M, Cao E, Julius D, Cheng Y. Structure of the TRPV1 ion channel determined by electron cryo-microscopy. Nature. 2013;504: 107–112. pmid:24305160
  7. 7. Opella SJ. Structure determination of membrane proteins in their native phospholipid bilayer environment by rotationally aligned solid-state NMR spectroscopy. Acc Chem Res. 2013;46: 2145–2153. pmid:23829871
  8. 8. Opella SJ, Marassi FM. Structure determination of membrane proteins by NMR spectroscopy. Chem Rev. 2004;104: 3587–3606. pmid:15303829
  9. 9. Vasquez V, Sotomayor M, Cordero-Morales J, Schulten K, Perozo E. A structural mechanism for MscS gating in lipid bilayers. Science. 2008;321: 1210–1214. pmid:18755978
  10. 10. Chakrapani S, Cuello LG, Cortes DM, Perozo E. Structural dynamics of an isolated voltage-sensor domain in a lipid bilayer. Structure. 2008;16: 398–409. pmid:18334215
  11. 11. Chakrapani S, Sompornpisut P, Intharathep P, Roux B, Perozo E. The activated state of a sodium channel voltage sensor in a membrane environment. Proc Natl Acad Sci USA. 2010;107: 5435–5440. pmid:20207950
  12. 12. Li Q, Wanderling S, Sompornpisut P, Perozo E. Structural basis of lipid-driven conformational transitions in the KvAP voltage-sensing domain. Nat Struct Mol Biol. 2014;21: 160–166. pmid:24413055
  13. 13. Borbat PP, McHaourab HS, Freed JH. Protein structure determination using long-distance constraints from double-quantum coherence ESR: study of T4 lysozyme. J Am Chem Soc. 2002;124: 5304–5314. pmid:11996571
  14. 14. Borbat PP, Freed JH. Measuring distances by pulsed dipolar ESR spectroscopy: spin-labeled histidine kinases. Methods Enzymol. 2007;423: 52–116. pmid:17609127
  15. 15. Jones DH, Barber KR, VanDerLoo EW, Grant CW. Epidermal growth factor receptor transmembrane domain: 2H NMR implications for orientation and motion in a bilayer environment. Biochemistry. 1998;37: 16780–16787. pmid:9843449
  16. 16. van der Wel PC, Strandberg E, Killian JA, Koeppe RE 2nd. Geometry and intrinsic tilt of a tryptophan-anchored transmembrane alpha-helix determined by (2)H NMR. Biophys J. 2002;83: 1479–1488. pmid:12202373
  17. 17. Alexander N, Bortolus M, Al-Mestarihi A, McHaourab H, Meiler J. De novo high-resolution protein structure determination from sparse spin-labeling EPR data. Structure. 2008;16: 181–195. pmid:18275810
  18. 18. Im W, Jo S, Kim T. An ensemble dynamics approach to decipher solid-state NMR observables of membrane proteins. Biochim Biophys Acta. 2012;1818: 252–262. pmid:21851810
  19. 19. Bhatnagar J, Freed JH, Crane BR. Rigid body refinement of protein complexes with long-range distance restraints from pulsed dipolar ESR. Methods Enzymol. 2007;423: 117–133. pmid:17609128
  20. 20. Sompornpisut P, Roux B, Perozo E. Structural refinement of membrane proteins by restrained molecular dynamics and solvent accessibility data. Biophys J. 2008;95: 5349–5361. pmid:18676641
  21. 21. Jo S, Im W. Transmembrane helix orientation and dynamics: insights from ensemble dynamics with solid-state NMR observables. Biophys J. 2011;100: 2913–2921. pmid:21689524
  22. 22. Pathak MM, Yarov-Yarovoy V, Agarwal G, Roux B, Barth P, Kohout S, et al. Closing in on the resting state of the Shaker K+ channel. Neuron. 2007;56: 124–140. pmid:17920020
  23. 23. Vargas E, Bezanilla F, Roux B. In search of a consensus model of the resting state of a voltage-sensing domain. Neuron. 2011;72: 713–720. pmid:22153369
  24. 24. Rulisek L, Vondrasek J. Coordination geometries of selected transition metal ions (Co2+, Ni2+, Cu2+, Zn2+, Cd2+, and Hg2+) in metalloproteins. J Inorg Biochem. 1998;71: 115–127. pmid:9833317
  25. 25. Harding MM. Small revisions to predicted distances around metal sites in proteins. Acta Crystallogr D Biol Crystallogr. 2006;62: 678–682. pmid:16699196
  26. 26. del Camino D, Holmgren M, Liu Y, Yellen G. Blocker protection in the pore of a voltage-gated K+ channel and its structural implications. Nature. 2000;403: 321–325. pmid:10659852
  27. 27. del Camino D, Yellen G. Tight steric closure at the intracellular activation gate of a voltage-gated K+ channel. Neuron. 2001;32: 649–656. pmid:11719205
  28. 28. Campos FV, Chanda B, Roux B, Bezanilla F. Two atomic constraints unambiguously position the S4 segment relative to S1 and S2 segments in the closed state of Shaker K channel. Proc Natl Acad Sci USA. 2007;104: 7904–7909. pmid:17470814
  29. 29. Lin MC, Abramson J, Papazian DM. Transfer of ion binding site from ether-a-go-go to Shaker: Mg2+ binds to resting state to modulate channel opening. J Gen Physiol. 2010;135: 415–431. pmid:20385745
  30. 30. Tao X, Lee A, Limapichat W, Dougherty DA, MacKinnon R. A gating charge transfer center in voltage sensors. Science. 2010;328: 67–73. pmid:20360102
  31. 31. Lin MC, Hsieh JY, Mock AF, Papazian DM. R1 in the Shaker S4 occupies the gating charge transfer center in the resting state. J Gen Physiol. 2011;138: 155–163. pmid:21788609
  32. 32. Henrion U, Renhorn J, Borjesson SI, Nelson EM, Schwaiger CS, Bjelkmar P, et al. Tracking a complete voltage-sensor cycle with metal-ion bridges. Proc Natl Acad Sci USA. 2012;109: 8552–8557. pmid:22538811
  33. 33. Bezanilla F. How membrane proteins sense voltage. Nat Rev Mol Cell Biol. 2008;9: 323–332. pmid:18354422
  34. 34. Yarov-Yarovoy V, Baker D, Catterall WA. Voltage sensor conformations in the open and closed states in ROSETTA structural models of K+ channels. Proc Natl Acad Sci USA. 2006;103: 7292–7297. pmid:16648251
  35. 35. Catterall WA. Molecular properties of voltage-sensitive sodium channels. Annu Rev Biochem. 1986;55: 953–985. pmid:2427018
  36. 36. Catterall WA. Voltage-dependent gating of sodium channels: correlating structure and function. Trends Neurosci. 1986;9: 7–10.
  37. 37. Guy HR, Seetharamulu P. Molecular model of the action potential sodium channel. Proc Natl Acad Sci USA. 1986;83: 508–512. pmid:2417247
  38. 38. Long SB, Campbell EB, Mackinnon R. Crystal structure of a mammalian voltage-dependent Shaker family K+ channel. Science. 2005;309: 897–903. pmid:16002581
  39. 39. Long SB, Tao X, Campbell EB, MacKinnon R. Atomic structure of a voltage-dependent K+ channel in a lipid membrane-like environment. Nature. 2007;450: 376–382. pmid:18004376
  40. 40. Payandeh J, Scheuer T, Zheng N, Catterall WA. The crystal structure of a voltage-gated sodium channel. Nature. 2011;475: 353–358. pmid:21743477
  41. 41. Zhang X, Ren W, DeCaen P, Yan C, Tao X, Tang L, et al. Crystal structure of an orthologue of the NaChBac voltage-gated sodium channel. Nature. 2012;486: 130–134. pmid:22678295
  42. 42. Payandeh J, Gamal El-Din TM, Scheuer T, Zheng N, Catterall WA. Crystal structure of a voltage-gated sodium channel in two potentially inactivated states. Nature. 2012;486: 135–139. pmid:22678296
  43. 43. Vargas E, Yarov-Yarovoy V, Khalili-Araghi F, Catterall WA, Klein ML, Tarek M, et al. An emerging consensus on voltage-dependent gating from computational modeling and molecular dynamics simulations. J Gen Physiol. 2012;140: 587–594. pmid:23183694
  44. 44. Khalili-Araghi F, Jogini V, Yarov-Yarovoy V, Tajkhorshid E, Roux B, Schulten K. Calculation of the gating charge for the Kv1.2 voltage-activated potassium channel. Biophys J. 2010;98: 2189–2198. pmid:20483327
  45. 45. Jensen MO, Jogini V, Borhani DW, Leffler AE, Dror RO, Shaw DE. Mechanism of voltage gating in potassium channels. Science. 2012;336: 229–233. pmid:22499946
  46. 46. Yarov-Yarovoy V, DeCaen PG, Westenbroek RE, Pan CY, Scheuer T, Baker D, et al. Structural basis for gating charge movement in the voltage sensor of a sodium channel. Proc Natl Acad Sci USA. 2012;109: E93–E102. pmid:22160714
  47. 47. Phillips JC, Braun R, Wang W, Gumbart J, Tajkhorshid E, Villa E, et al. Scalable molecular dynamics with NAMD. J Comput Chem. 2005;26: 1781–1802. pmid:16222654
  48. 48. Humphrey W, Dalke A, Schulten K. VMD: visual molecular dynamics. J Mol Graphics. 1996;14: 33–38.
  49. 49. Roux B, Islam SM. Restrained-ensemble molecular dynamics simulations based on distance histograms from double electron-electron resonance spectroscopy. J Phys Chem B. 2013;117: 4733–4739. pmid:23510121
  50. 50. Islam SM, Stein RA, McHaourab HS, Roux B. Structural refinement from restrained-ensemble simulations based on EPR/DEER data: application to T4 lysozyme. J Phys Chem B. 2013;117: 4740–4754. pmid:23510103
  51. 51. Islam SM, Roux B. Simulating the distance distribution between spin-labels attached to proteins. J Phys Chem B. 2015;119: 3901–3911. pmid:25645890
  52. 52. Mangani S, Carloni P, Orioli P. Crystal structure of the complex between carboxypeptidase A and the biproduct analog inhibitor L-benzylsuccinate at 2.0 Å resolution. J Mol Biol. 1992;223: 573–578. pmid:1738164
  53. 53. Kumaran D, Eswaramoorthy S, Furey W, Sax M, Swaminathan S. Structure of staphylococcal enterotoxin C2 at various pH levels. Acta Crystallogr D Biol Crystallogr. 2001;57: 1270–1275. pmid:11526318
  54. 54. Papageorgiou AC, Acharya KR, Shapiro R, Passalacqua EF, Brehm RD, Tranter HS. Crystal structure of the superantigen enterotoxin C2 from Staphylococcus aureus reveals a zinc-binding site. Structure. 1995;3: 769–779. pmid:7582894
  55. 55. Christianson DW, Lipscomb WN. X-ray crystallographic investigation of substrate binding to carboxypeptidase A at subzero temperature. Proc Natl Acad Sci USA. 1986;83: 7568–7572. pmid:3463986
  56. 56. Russell RB, Barton GJ. Multiple protein sequence alignment from tertiary structure comparison: assignment of global and residue confidence levels. Proteins. 1992;14: 309–323. pmid:1409577
  57. 57. Essmann U, Perera L, Berkowitz ML, Darden T, Lee H, Pedersen LG. A smooth particle mesh ewald method. J Chem Phys. 1995;103: 8577–8593.
  58. 58. Huang J, MacKerell AD Jr. CHARMM36 all-atom additive protein force field: validation based on comparison to NMR data. J Comput Chem. 2013;34: 2135–2145. pmid:23832629
  59. 59. MacKerell AD, Bashford D, Bellott M, Dunbrack RL, Evanseck JD, Field MJ, et al. All-atom empirical potential for molecular modeling and dynamics studies of proteins. J Phys Chem B. 1998;102: 3586–3616. pmid:24889800
  60. 60. MacKerell AD Jr, Feig M, Brooks CL 3rd. Improved treatment of the protein backbone in empirical force fields. J Am Chem Soc. 2004;126: 698–699. pmid:14733527
  61. 61. Jorgensen WL, Chandrasekhar J, Madura JD, Impey RW, Klein ML. Comparison of simple potential functions for simulating liquid water. J Chem Phys. 1983;79: 926–935.
  62. 62. Ryckaert JP, Ciccotti G, Berendsen HJC. Numerical integration of cartesian equations of motion of a system with constraints: molecular dynamics of n-alkanes. J Comput Phys. 1977;23: 327–341.
  63. 63. Martyna GJ, Tobias DJ, Klein ML. Constant pressure molecular dynamics algorithms. J Chem Phys. 1994;101: 4177–4189.
  64. 64. Feller SE, Zhang YH, Pastor RW, Brooks BR. Constant pressure molecular dynamics simulation: the Langevin piston method. J Chem Phys. 1995;103: 4613–4621.