PMCCPMCCPMCC

Search tips
Search criteria 

Advanced

 
Logo of nihpaAbout Author manuscriptsSubmit a manuscriptHHS Public Access; Author Manuscript; Accepted for publication in peer reviewed journal;
 
Biochemistry. Author manuscript; available in PMC 2010 December 22.
Published in final edited form as:
PMCID: PMC2795120
NIHMSID: NIHMS159029

Apo and Nickel-bound Forms of the Pyrococcus Horikoshii Species of the Metalloregulatory Protein: NikR Characterized by Molecular Dynamics Simulations

Abstract

NikR is a homotetrameric nickel regulatory protein whose binding to free Ni2+ increases its binding affinity for a gene that codes for a nickel transporter protein. It is comprised of a tetrameric nickel-binding domain, flanked by two dimeric DNA-binding domains. Though x-ray crystallography data for various species (Escherichia coli, Heliobacter pylori, and Pyrococcus horikoshii) of NikR reveal large conformational differences between Nickel-bound, DNA-bound, and unbound forms, transitions between them have never been observed. We have run all-atom molecular dynamics simulations of three forms of the Pyrococcus horikoshii species of NikR including two apo-forms and one nickel-bound form. Though all 552 residues of this species occur naturally, quantum-mechanics-based force-field parameterization was required to accurately represent the four nickel-centers in the nickel-bound form. Global conformational analysis of the three 100-ns-long simulations indicates slow conformational kinetics and independent DNA binding domain motion. Correlation and flexibility analysis revealed regions of high structural and dynamical importance. A striking relationship was observed between regions with high levels of structural importance and regions with known biological importance. Mutation of key regions of P. horikoshii and analogous regions in both E. coli and H. pylori are suggested that might inhibit DNA-binding activity while not affecting nickel-binding.

Keywords: molecular dynamics, nickel regulation, protein allostery, metallo-protein, NikR

The understanding of metal ion homeostasis at the molecular-level is critical to our understanding of how toxic metal ions are delivered safely and specifically to the appropriate location within cells.13 There is a significant structural, biochemical and medical literature on the homeostasis of Cu, for example, due to its role in several neurodegenerative disorders (Alzheimer’s and amyotrophic lateral sclerosis (Lou Gehrig’ disease), for example).1, 2, 48 Despite its inherent toxicity, many microorganisms require cellular nickel as an enzyme co-factor.9, 10 Escherichia coli and Heliobacter pylori, for example, utilize nickel for [NiFe]-hydrogenases and urease, which are important for microbial energy metabolism and creation of ammonia as nitrogen source, respectively.10, 11 NikR is a DNA-binding metalloregulatory protein that utilizes nickel-responsive allostery to modulate cellular nickel concentrations in various ways.3, 12, 13 In E. coli, NikR detects cellular nickel levels and selectively binds and represses the operon which encodes NikABCDE, an ABC-type membrane transporter12, 14 effectively reducing nickel import. In H. pylori, NikR, activated by high nickel levels, represses expression of nickel-dependent proteins, such as the HoxN type permease (NixA) which imports nickel, while stimulating activation of other proteins.15 Hence, NikR senses intracellular levels of Ni and shuts downs the production of the transporter when the level of Ni is appropriate for the cell to function.12

Allosteric mechanism

It has been found that NikR’s affinity for DNA is dependent on nickel concentration.16, 17 NikR contains both low (30 nM) and high-affinity (20 pM) nickel-binding sites; when the high (low) affinity nickel sites become occupied, DNA affinity increases to nanomolar (picomolar) levels.17 Carrington et al. found that nickel binds to the high affinity sites in a unique square-planar geometry, but changes to 6-coordinate when bound to DNA.18 Several species of NikR have been crystallized in various forms by Drennan and co-workers (E. coli NikR – EcNikR)19, 20 Terradot and co-workers (H. pylori NikR - HpNikR)21 and Chivers and Tahirov (Pyrococcus horikoshii NikR - PhNikR).22 The overall structure consistently shows NikR to have a central tetrameric metal-binding domain (MBD) constructed of a βαββαβ ferredoxin-like fold binding nickel in square-planar geometry flanked by two dimeric ribbon-helix-helix (βαα) DNA binding domains (DBDs) (see Figure 1).

Figure 1
A rendering of PhNikR. The DNA binding domains (βαα) surround the central metal binding domain (βαββαβ). Metallic spheres represent bound nickels.

There are, however, some unique characteristics between different species’ crystal structures. The relatively stable MBD of NikR and two rotatable DBDs of NikR can adapt roughly three different main conformations (straight, trans or cis) whose stability seems to be related to both nickel and DNA binding. We will refer to the different forms by the presence of nickel, A (apo) for nickel-free, H (holo) for nickel-bound, and the global conformation, T for trans-like, C for cis-like, and S for straight. Figure 2 shows sketches of various forms of NikR. Schreiter et al. crystallized EcNikR in AS, HS and HC (DNA-bound) forms.19, 20, 23 They found that only the nickel-bound forms had a well-ordered α3 helix (a helix containing one of the four nickel-binding residues). This lead to the proposition that the stabilization of the α3 helix assists in DNA binding, a hypothesis supported by comparative metal binding studies.24 Chivers and Tahirov crystallized five forms of P. horikoshii including AS, HS, and three forms of HT: one with the four high affinity nickel-binding sites occupied and PhNikRs unique additional four auxiliary high affinity sites bound, one with all high and low sites bound, and one with all sites bound in the presence of free phosphate ions (to simulate DNA).22 The AS form, similar to EcNikR, was found to have a highly mobile and disordered α3-helix. Interestingly, the form with phosphate ions is oriented in a trans-form that should not fit into the corresponding DNA operon assuming that the DNA binding is similar to that of EcNikR which binds in a cis-form. As Chivers and Tahirov noted, this form can be converted into a cis-form if one half (any two chains that share a DBD) is rotated by 180 degrees about the center axis (turning it from a trans to a cis).22 The crystal forms found by Dian et al. of HpNikR were the most unique set of the three species crystallized.21 They found a unique AT structure, which roughly resembles the HT form of PhNikR but lacks dimeric symmetry. Also, HpNikR was never found with all four high-affinity nickel-binding sites occupied, but rather with a single high-affinity nickel site bound per dimer, an intermediate site, or an external site. It is possible that this peculiarity could have been caused by the low pH conditions of the crystallization,21 as noted by Zambelli et al.25

Figure 2
Sketch of various possible forms of NikR. Rectangles represent metal binding domains (each color represents a different chain, purple and lime represent intertwined chains). Trapezoids represent the DNA binding domains (each has contributions from two ...

These various permutations of conformations and nickel and DNA bound states both provide insight and raise questions about NikR’s allosteric mechanism. What role does nickel binding have in the structural preference of NikR? Does nickel-binding preferentially stabilize one conformation over another otherwise in preexisting equilibrium?26 What residues are important in promoting these conformational changes?

To elucidate the role of nickel binding on conformational states and dynamics of NikR we simulated various forms using all-atom molecular dynamics (MD) simulations. MD simulations of biological molecules have often been utilized to elucidate biological mechanisms.27 Few computational studies of any species of NikR have been reported. Cui and Merz utilized a course-grained Gaussian-network-model to find natural global fluctuations in PhNikR.28 They found that natural fluctuations lead to transitions between the major conformational states found in NikR crystal structures (between straight, trans, and cis). This coarse-grained study supports the pre-existing equilibrium hypothesis – that all the available conformational states are available to each bound state.26 Bradley et al. utilized29 an 80 ns all-atom explicit solvent simulation of EcNikR starting from the AT state (PDB code 1Q5V). In that work, a potential interaction network between the metal-binding and DNA binding residues was highlighted using a novel correlation-matrix-based residue clustering technique.29 However, the metal binding states we not examined. During the preparation of this manuscript, a computational study by Phillips et al. was published regarding NikR metal binding specificity based on free energy simulations.30 This study utilized relatively stable, DNA-bound states of EcNikR for Poisson-Boltzmann Thermodynamic Integration calculation of free energy differences between different bound ions.

Herein, we examine the dynamics of three 100 ns simulations beginning from three forms of PhNikR in solution. The Ph species was chosen because the crystal structures available had the least amount of missing atoms. We chose the AS, HT (with only four nickels bound) and AT forms (PDB codes: 2BJ3, 2BJ7 (with 4 primary nickel ions), 2BJ7 (without nickel ions)). Our goal is to highlight the features of the dynamics of each form to compare and contrast any apparent differences between them. The analyses reported herein include: RMSD analysis to look at conformational changes of various regions of interest, center-of-mass angular analysis to measure global conformational change, correlation analysis to examine inter-residue communication, and flexibility analysis to identify residues necessary for hinge-motion.

Materials and Methods

Simulation details

We chose to simulate three forms of PhNikR22: the apo-straight form, taken from PDB ID 2BJ3 (PhNikR-AS), the holo-trans form with only four nickels simulated, taken from PDB ID 2BJ7 (PhNikR-HT), and an apo-trans from where all nickels were stripped from the original pdb structure 2BJ7 (PhNikR-AT). An implicit generalized Born solvation model (GB(OBC))31 was chosen to reduce the system size as well as the solvent viscosity to allow for faster kinetics.32 Recently, Hornak et al. successfully utilized implicit solvation to witness flap-motion in a simulation of HIV-1 protease.33,34 In our simulations, each form was taken from the appropriate PDB structure with missing residues built by the TLEAP program, which is part of the AMBER suite.35 The PhNikR-HT form required additional modeling in order to properly represent the nickel centers (see the “Force Field Parameters” section). The proteins were simulated in AMBER using the ff99SB parameter set.36 Hydrogen vibrations were restrained using the SHAKE algorithm allowing a 2 fs timestep. Each form was minimized in AMBER10 using the SANDER module, then relaxed at 100K for 10 ps. The production runs were 100 ns long at 300K for each form using the PMEMD module. Temperature was regulated by utilization of Langevin dynamics using a 1.0 ps−1 solvent coupling constant. Successive restart runs utilized different pseudorandom seeds as recommended by Sindhikara et al.37 Analysis of the trajectories was performed using the PTRAJ module.

Force field parameters

All 138×4 of the amino acid residues in Ph apo-NikR occur naturally and thus force field parameters are available in various parameter sets (ff99SB in our case36). For the apo forms, all histidine protonation states were chosen to agree with the most likely protonation state at the experimental pH of 7.5 and are thus delta-protonated. Since in the apo form no histidines are hydrogen bound any other residues, the difference between the delta and epsilon (the second most likely state) protonation will have a limited impact, if at all, to global dynamics.

For the nickel-bound form, there are several approaches to metal modeling that we could have chosen. These approaches can roughly be classified into three categories: quantum mechanical, non-bonded, and bonded. Quantum mechanical treatment or a mixed QM/MM treatment of metals has the most potential for accuracy, but its huge computational cost deters usage for large systems. A non-bonded approach would treat the nickel as a free ion whose interaction with the protein is strictly through van der Waals and Coulomb forces. This is the approach Phillips et al. in their metal preference study.30 This approach is delicate to parameterize, does not enforce a particular geometry, and may ignore any covalent bond characteristics.3842 Finally, a bonded approach utilizes all terms in the typical all-atom force field to treat the metal ion. Parameterization of this approach can be done by combination of experimental and quantum mechanical methods. This method is the most restrictive geometrically, but is useful when modeling stable metal sites. We chose the bonded approach since NikR is too large to treat quantum mechanically, and we know precisely what the Ni coordination state is in NikR.

New force field parameters were created and used for both the bound nickel ions, and the ligating residues. The MTK++ package was used to expedite the force field parameterization process.43 First, a single copy of one of the nickels plus the four ligating sidechains were isolated from the holo structure (pdb code: 2BJ7) then capped with CH3 groups (in the position of the original ). Figures 3a and b show representations of this “small” metal-cluster model.

Figure 3
Representations of the “small” metal cluster. Atoms included in this model were used to obtain force constants via DFT calculations. a) A ball-and-stick representation. Labels indicate residue identity and chain ID b) A 2D representation. ...

The histidine protonation states were taken from experiment: H91-A and H78-D are epsilon-protonated, and H89-A is delta-protonated with the remaining sites occupied by the Nickel ion. The resultant cluster was then minimized using the TPSSTPSS/TZVP level of DFT theory with Gaussian 03.44 This method and basis has been shown to model nickel compounds well.45 The minimization did not alter the square-planar geometry (though use of smaller basis sets did, data not shown). Following the minimization, force constant calculations were performed using the same methods.

To obtain the partial charges on all the atoms, a larger cluster was selected from the pdb file, this time including backbone atoms, but capped with Acetyl and N-Methyl groups in appropriate positions. Figures 4a and b show representations of this “large” cluster.

Figure 4
Representations of the “large” metal cluster. Atoms included in this model were used to obtain partial atomic charges via DFT. a) A ball-and-stick representation. b) A 2D representation.

Charges were obtained using the RESP algorithm4653 on this model cluster. van der Waals parameters were not adjusted for the ligating residues. For the Nickel, thermodynamic-integration simulations were performed for Ni2+ + 2Cl- + nH2O → nH2O for various ε and r values for nickel. The combination of r and ε that most closely resulted in the experimental free energy of hydration was used for nickel.54 The new set of parameters for the nickel ions and ligating residues were then mapped to the four primary nickel-binding sites in the holo form of NikR. The ff99SB parameter set36 was used for all other residues. Note that although no restraints were used to maintain the planarity of the nickel center the planarity was maintained throughout the course of the simulations described below (data not shown).

Table 1 shows the partial charges of the nickel center along with those of the nickel-bound atoms. For the nickel ion, we used a van der Waals radius of 1.10Å and a well-depth of 0.013 kcal/mol in all instances. All additional parameters are available in the supporting information.

Table 1
Classification of structure based on angles between centers of mass of different domains.

Results and Discussion

Simulation stability

To evaluate the energetic relaxation of the simulations, the potential energy of each system was examined. Figure 5 shows a block-averaged time-series of the potential energy. AS refers to the straight nickel-free structure (pdb ID 2BJ3), AT the trans nickel-free structure (pdb ID 2BJ7), and HT is the nickel-bound trans structure (pdb ID 2BJ7).

Figure 5
Averaged potential energy time series for three simulations of NikR. The smoothing window size is 50 ps.

The potential energy for all simulations had standard deviations below 80 kcal/mol for the last 90% of the simulation. This is consistent with the fluctuation in the simulation by Bradley et al.29 whose potential energy fluctuated by about 200 kcal/mol. It should be noted, however, that these systems potential landscapes are quite different considering the different solvent representations. The potential energy drops off slightly for the first 10 ns. After that, all three simulations fluctuate about each of their respective averages. Thus for the final 90 ns, each simulation is likely in thermal equilibrium.

Global conformational dynamics

Observing the global RMSD vs. time can give a first glimpse as to whether major conformational changes occurred during the simulation. A relatively static RMSD compared to the initial structure suggests a stable conformation while a dynamic RMSD suggests significant motion (thought it gives little insight to the specifics of that motion). Figures 6a–b show the RMSD of Cα atoms for key subsets of the NikR structure. Figure 6a shows the Cα RMSD against the initial structure versus time for all three simulations (AS, AT, and HT). All three simulations deviate more than 2Å within the first nanosecond. Though none of the simulations are stable over the entire 100 ns of simulation, they each maintain a steady RMSD between transitions. Both of the simulations whose initial structures came from forms found via crystallography (AS and HT) maintain states for nearly 60 ns. The AT simulation only stays steady for ~40 ns and experiences the highest RMSD change. The magnitude of fluctuations is similar to that of the 80 ns explicit solvent E. Coli simulation of Bradley et al.29 whose backbone RMSD, after relaxation, fluctuated between 2.5 and 5Å. The greater fluctuation observed in the simulations herein is likely due to the lack of explicit solvent molecules in the GB simulation. This inhomogeneous behavior suggests that none of the three systems have attained conformational equilibrium within the limit of the 100 ns timescale used in our simulations. Thus, we conclude that our analyses of the resultant trajectories only provide qualitative insights into these fascinating systems.

Figure 6
RMSD plots for Cα atoms of selected regions of NikR. a) Cα RMSD for entire protein in the AS, AT and HS simulations b) Cα RMSD for the metal binding domain for all three simulations.

The x-ray structures available for NikR suggest that the central MBD and two flanking DBDs may act similarly to rigid bodies. If each of these domains were rigid, the Cα RMSD would be relatively low and stable throughout the simulation. Figure 4b shows the Cα RMSD against the first frame for only the MBD of the protein (Residues 57 to 131 of each chain). Both the AS and HT simulations jump initially and stay in a steady state between 1.5 and 2Å from the initial structure. Thus, for these two forms, the central MBD is relatively rigid throughout the simulation. For the AS and HT forms, the large conformational change around 60 ns indicated by the full protein RMSD in Figure 6a must not be directly due to changes within the MBD. The AT simulation jumps initially to 3Å then rises to 4Å between the 40th and 60th ns. This larger jump may be due to the instability of the initial nickel-free AT structure and may be the source of the global conformational change indicated in Figure 6a. Since the MBD is relatively stable for both the AS and HT simulations, the conformational changes must come from the DBD motion. RMSD plots for DBDs can be found in the SI.

To observe the rigid body motion of the DBDs relative to the MBD, motion was observed in the context of the transition between the straight, trans, and cis forms referred to in Figure 2. Consider the centers of mass (COMs) of the four major domains: DBD from chains A and B (DBDAB), MBD from AB (MBDAB), MBD from CD (MBDCD), and DBD from CD (DBDCD). We can classify the shape of the protein based on these four points as in Table 2.

Table 2
Partial charges for the nickel center along with each ligating atom. Partial charges for each atom in the ligating residues are provided in the supporting information.

Figures 7a and b show the angles between both (DBDAB, MBDAB, and MBDCD) and (MBDAB, MBDCD, and DBDCD) for all three forms. The AS form, clearly maintains angles near 180 degrees, indicating that the COMs of the four domains are roughly linear. Conversely, the AT and HT forms maintain an angle of 150 degrees or less -- indicating a bend in the structure. While fluctuations do exist, no major conformational changes are evident from this analysis other than the HT CD DBD, which drifts near 160 degrees. This could indicate that this DBD is “straightening out”. Interestingly, the flapping motion seems to be independent, suggesting that the global conformational rearrangement may act in a step-wise fashion, rather than in concert.

Figure 7
Angles between centers of mass of adjacent domains vs time. a) Angle made by DBDAB, MBDAB, and MBDCD. b) Angle made by MBDAB, MBDCD, and DBDCD. c) Dihedral angle between DBDAB, MBDAB, MBDCD, and DBDCD.

Figure 7c shows the dihedral angle between these four points for each simulation. Unsurprisingly, the AS domain-COM dihedral angle fluctuates wildly across all angles since it is nearly linear and thus small motions can cause large angular changes. Conversely, the AT and HT, both starting at less than − 150 degrees, slowly drift up throughout the simulation. This could be simply relaxation, or more significantly, a slow conformational change to the cis form.

These Cα RMSD and domain COM angular analyses indicate significant conformational motion, but no complete transitions between the major states observed in the crystallographic experiments. However, the CD DBD of the HT simulation begins to straighten out within the timescale of the simulation. This motion suggests that conformational transitions may be observable within an order of magnitude longer simulation (~1 µs).

Features of key structural segments

Though the relative motion of the large domains reflects global conformational changes, the behavior of the segments that comprise and connect these regions can reflect the communication between them. We identified structural segments with key functionalities by analyzing both residue-level distance correlation and flexibility.

(A) Residue Correlation

To elucidate inter-residue dynamical relationships, distance correlation analysis was performed on all Cα atoms for each simulation with the PTRAJ module from AMBER10 (see Equation 1)

C(i,j)=Δxi(tk)Δxi(tk)Δxi(tk)21/2Δxj(tk)21/2
(1)

The correlation matrices for each simulation (found in SM) contain much information but are difficult to parse visually. Bradley et al. reported visually similar correlation plots for their EcNikR simulation. In order to quantify the relevance of specific residues, we summed the absolute value of the correlation matrix along either a row or column for each residue summed over the four chains as in Equation 2

cRELi=chainsj=1M|ci,j|
(2)

Here M is the number of residues per chain. Figure 8a shows this correlation relevance summed over chains (cREL). The colored area just above the x-axis indicates biological information including binding residues and highly conserved residues across 12 selected species of NikR.5565 We chose to use the same 12 species examined in the original PhNikR crystal structure paper by Chivers et al.22 The DNA binding residues include those involved in nonspecific contacts. The graphs suggest that each form exhibits similar behavior, residue-to-residue, and that the residues of AS and HT forms are much more coupled. Interestingly, for most of the DBD (approximately first 38 residues) for these two forms, the chain-independent correlation relevance values are nearly identical. In the MBD, the AS is more correlated than the HT form (and even more so than the AT). Additionally, the graph suggests that residues with known biological function tend to have higher correlation relevance with the exception of Gly-67, which may require more flexibility.

Figure 8
Correlation relevance summed over chains (cREL) for all three simulations. a) cREL vs residue number. Residue identities are labeled at peak values. Bottom bars indicate residues with high conservation across species, or with either nickel or DNA binding ...

Figure 8b shows the values cREL of the HT simulation (though all simulations are qualitatively similar in this respect) applied to a color scheme (red is the highest, blue is the lowest). Remarkably, both the segments of secondary structure involved in nickel-binding (α3 and β3) and those involved in DNA-binding (β1) have the highest value of cREL. This suggests that these regions are the most important structural regions since they coordinate motion the strongest with other residues. Mutation of the residues that would disrupt the structural stability of these segments would likely inhibit DNA-binding. Previously performed mutation studies of some of the residues for various species in these segments have already shown the importance of these segments in DNA and nickel-binding.13, 18, 21 However, the DNA binding characteristics of each species is unique; thus, mutation studies should be performed on different but analogous residues. Additional correlation analysis can be found in the SM. We note that Bradley et al. utilized a different correlation-based measure to determine interaction groups rather than individual residue relevance as shown here.29

(B) Residue Flexibility

Flexibility can be a necessary attribute in key regions of a protein to allow it to not only fold but also undergo functional conformational changes. For example, flexibility in hinge regions is advantageous.6668 To identify regions of high flexibility in NikR, each residue was characterized by the standard deviation of the distance between the Cα atoms of its nearest neighbors summed over chains (see Equation 3)

cFLXi=chainsσ(|xi1xi+1|)
(3)

Here [x with verticle line above]i represents the position of the Cα atom of the ith residue. Figure 9a shows cFLX for each simulation. The flexibility seems to be almost identical for all three simulations suggesting the form has little effect on the residue-level motion. High flexibility is seen at the two termini as well as in three regions centered on Glu-28, Gly-50, and Asp-68 as labeled. It is unsurprising that the correlation relevance, cREL, and flexibility, cFLX, do not correlate with each other. cREL corresponds more to structural stability while cFLX corresponds to instability. We note that there is no correlation between the X-ray B-factor and this measure of flexibility (data not shown). This loop flexibility may play a key role in allowing the allosteric conformational change required for NikR’s function. A mutation to a residue within a region of high flexibility to increase stability should decrease the biological function.

Figure 9
Residue flexibility, cFLX, for all three simulations. a) cFLX vs. residue number. Residue identities are labeled at peak values. Bottom bars indicate residues with high conservation across species, or with either nickel or DNA binding function as labeled. ...

Figure 9b shows the values cFLX of the HT simulation in the same manner as in 8b applied to a color scheme (red is the highest, blue is the lowest). The figure shows that the highest flexibility regions are in only a few of the loop regions, particularly the α1-α2, α2-β2, and α3-β3 loops (these correspond to the three peaks in Figure 9a. The analogous flexible loops for E. coli are Gly-24 to Asn-27, Gln-47 to Gln-51 and Glu-70 to Asp-73, and for H. pylori are Gly-33 to Ser-36, Glu-56 to Asn-60, and His-72 to Glu-85. While all three of these loop regions may have a distinct biological function, the α2-β2 loop is the most interesting in the context of the global conformational changes NikR undergoes since it lies at the interface between the DBD and MBD regions. The hinge-motion of this loop is key for the switching between the cis, trans, and straight conformations. Interestingly, the sequence conservation in this region is low.5565 This may indicate that these residues are necessary in DNA-binding specificity.

Since the α2-β2 loop is a hinge that allows the global conformational change required for NikR to bind to DNA (converting to the cis-form), and in the wild type it is highly flexible, mutations to this region that restrict motion should inhibit DNA binding activity. Two key features of this loop are highlighted in Figure 10. The first is the Gly-50 residue that most known species of NikR contain (it is absent in the α2-β2 loop of H. pylori).5565 The mutation of this glycine to a conformationally less flexible residue such as a proline might reduce NikRs ability to change conformation to accommodate DNA binding. Secondly, there is a hydrophobic interaction between Val-49 in the α2-β2 loop and Met-1 in the β5 strand. The reduction of the length of either of these sidechains could reduce the flexibility as well. For example, mutation of the Met-1 to a Val or Val-49 to Ala would keep the hydrophobic interaction but would tie-down the loop to restrict motion. Similar mutations on other species of NikR would likely produce the same effect. For example, in HpNikR, this interaction is coulombic rather than hydrophobic (Lys-31 with Glu-56) as the crystal structure (PDBID: 2CA9) suggests. For E. coli, the interactions is not obvious from the crystal structure (PDBID: 2HZA), however, the study by Bradley et al. did indeed recognize a non-covalent interaction in this general region (in that a unified cluster was recognized containing both sides of this region).

Figure 10
Close-up of α2-β2 loop. Key residues are represented by VDW volume spheres and labeled appropriately. Gly-50 provides flexibility. Shortening the sidechains while keeping the hydrophobic interaction can reduce the Cα-Cα ...

Conclusions

NikR regulates cellular nickel concentration by responding to high nickel concentrations by suppressing nickel transport protein expression. Experimental data gives clues to its mechanism but without conclusive detail. We have performed three simulations to help elucidate the dynamics of NikR. The modeling of this system required special parameterization to accurately represent the nickel and its ligating residues. Our analyses of these simulations provide additional insight into the dynamical preference that different initial conformations and nickel-binding states have on the protein. While not observing global conformational changes on the computational time-scale available, several distinct features of NikR were observed. The flanking DBDs move independently of each other. Analysis of biologically relevant residues tend to influence motion across domains much more strongly than others residues.

The secondary structure segments directly pertaining to binding of either DNA or nickel were found to be structurally important by correlation analysis even in the absence of both of these ligands. Additionally, key flexible regions were highlighted especially the α2-β2 loop whose hinge-like motion allows functional global conformational change. Regions with either high correlation relevance or flexibility seem to be correlated with biologically important regions. We suggested some specific mutations to these residues for the P. horikoshii species and also identified analogous regions in E. coli and H. pylori species of NikR that we expect to have similar character.

With current computational capacity, simulating a global conformational change in a system like NikR is a significant challenge. Beyond brute force simulations, it is possible that advanced sampling algorithms could be applied to not only establish a meaningful reaction coordinate between NikR states, but also to calculate the free-energy change between states.

Supplementary Material

1_si_001

ACKNOWLEDGEMENTS

We thank Martin B. Peters for assistance with the metal modeling software (MTK++). We also thank the anonymous reviewer for excellent and thorough reviews and suggestions.

FUNDING INFORMATION

This work was supported by the National Institutes of Health (NIH) grant GM066689. Simulations were performed on the Teragrid “Kraken” cluster under grant #: MCA05S010 and the University of Florida HPC (UFHPC).

Footnotes

Supporting Information Available

Additional RMSD plots, correlation plots, fluctuation analysis, as well as force field parameters. This material is available free of charge via the Internet at http://pubs.acs.org

REFERENCES

1. Culotta VC, Yang M, O'Halloran TV. Activation of superoxide dismutases: Putting the metal to the pedal. Biochimica Et Biophysica Acta-Molecular Cell Research. 2006;1763(7):747–758. [PMC free article] [PubMed]
2. Finney LA, O'Halloran TV. Transition metal speciation in the cell: insights from the chemistry of metal ion receptors. Science. 2003;300(5621):931–937. [PubMed]
3. Dosanjh NS, Michel SLJ. Microbial nickel metalloregulation: NikRs for nickel ions. Curr Opin Chem Biol. 2006;10(2):123–130. [PubMed]
4. Valentine JS, Doucette PA, Potter SZ. Copper-zinc superoxide dismutase and amyotrophic lateral sclerosis. Annu Rev Biochem. 2005;74:563–593. [PubMed]
5. Rosenzweig AC. Copper delivery by metallochaperone proteins. Acc Chem Res. 2001;34(2):119–128. [PubMed]
6. Gaggelli E, Kozlowski H, Valensin D, Valensin G. Copper homeostasis and neurodegenerative disorders (Alzheimer's, prion, and Parkinson's diseases and amyotrophic lateral sclerosis) Chem Rev. 2006;106(6):1995–2044. [PubMed]
7. Beckman JS, Estevez AG, Crow JP, Barbeito L. Superoxide dismutase and the death of motoneurons in ALS. Trends Neurosci. 2001;24(11 Suppl):S15–S20. [PubMed]
8. Huffman DL, O'Halloran TV. Function, structure, and mechanism of intracellular copper trafficking proteins. Annu Rev Biochem. 2001;70:677–701. [PubMed]
9. Maroney MJ. Structure/function relationships in nickel metallobiochemistry. Curr Opin Chem Biol. 1999;3(2):188–199. [PubMed]
10. Mulrooney SB, Hausinger RP. Nickel uptake and utilization by microorganisms. FEMS Microbiol Rev. 2003;27(2–3):239–261. [PubMed]
11. Mobley HLT, Hausinger RP. Microbial Ureases - Significance, Regulation, and Molecular Characterization. Microbiological Reviews. 1989;53(1):85–108. [PMC free article] [PubMed]
12. De Pina K, Desjardin V, Mandrand-Berthelot MA, Giordano G, Wu LF. Isolation and characterization of the nikR gene encoding a nickel-responsive regulator in Escherichia coli. J Bacteriol. 1999;181(2):670–674. [PMC free article] [PubMed]
13. Chivers PT, Sauer RT. NikR is a ribbon-helix-helix DNA-binding protein. Protein Sci. 1999;8(11):2494–2500. [PubMed]
14. Chivers PT, Sauer RT. Regulation of high affinity nickel uptake in bacteria. Ni2+-Dependent interaction of NikR with wild-type and mutant operator sites. J Biol Chem. 2000;275(26):19735–19741. [PubMed]
15. Contreras M, Thiberge JM, Mandrand-Berthelot MA, Labigne A. Characterization of the roles of NikR, a nickel-responsive pleiotropic autoregulator of Helicobacter pylori. Mol Microbiol. 2003;49(4):947–963. [PubMed]
16. Chivers PT, Sauer RT. NikR repressor: High-affinity nickel binding to the C-terminal domain regulates binding to operator DNA. Chem Biol. 2002;9(10):1141–1148. [PubMed]
17. Bloom SL, Zamble DB. Metal-selective DNA-binding response of Escherichia coli NikR. Biochemistry. 2004;43(31):10029–10038. [PubMed]
18. Carrington PE, Chivers PT, Al-Mjeni F, Sauer RT, Maroney MJ. Nickel coordination is regulated by the DNA-bound state of NikR. Nature Structural Biology. 2003;10(2):126–130. [PubMed]
19. Schreiter ER, Sintchak MD, Guo YY, Chivers PT, Sauer RT, Drennan CL. Crystal structure of the nickel-responsive transcription factor NikR. Nature Structural Biology. 2003;10(10):794–799. [PubMed]
20. Schreiter ER, Wang SC, Zamble DB, Drennan CL. NikR-operator complex structure and the mechanism of repressor activation by metal ions. Proc Natl Acad Sci U S A. 2006;103(37):13676–13681. [PubMed]
21. Dian C, Schauer K, Kapp U, McSweeney SM, Labigne A, Terradot L. Structural basis of the nickel response in Helicobacter pylori: Crystal structures of HpNikR in Apo and nickel-bound states. J Mol Biol. 2006;361(4):715–730. [PubMed]
22. Chivers PT, Tahirov TH. Structure of Pyrococcus horikoshii NikR: Nickel sensing and implications for the regulation of DNA recognition. J Mol Biol. 2005;348(3):597–607. [PubMed]
23. Schreiter ER, Chivers PT, Drennan CL. Crystal structures of the Ni-dependent transcriptional repressor NikR from E. coli. FASEB J. 2004;18(8):C188–C188.
24. Phillips CM, Schreiter ER, Guo Y, Wang SC, Zamble DB, Drennan CL. Structural basis of the metal specificity for nickel regulatory protein NikR. Biochemistry. 2008;47(7):1938–1946. [PubMed]
25. Zambelli B, Bellucci M, Danielli A, Scarlato V, Ciurli S. The Ni2+ binding properties of Helicobacter pylori NikR. Chemical Communications. 2007;35:3649–3651. [PubMed]
26. Koshland DE, Nemethy G, Filmer D. Comparison of Experimental Binding Data and Theoretical Models in Proteins Containing Subunits. Biochemistry. 1966;5(1):365–385. [PubMed]
27. Karplus M, McCammon JA. Molecular dynamics simulations of biomolecules. Nature Structural Biology. 2002;9(9):646–652. [PubMed]
28. Cui G, Merz KM. The intrinsic dynamics and function of nickel-binding regulatory protein: Insights from elastic network analysis. Biophys J. 2008;94(10):3769–3778. [PubMed]
29. Bradley MJ, Chivers PT, Baker NA. Molecular dynamics simulation of the Escherichia coli NikR protein: Equilibrium conformational fluctuations reveal interdomain allosteric communication pathways. J Mol Biol. 2008;378(5):1155–1173. [PMC free article] [PubMed]
30. Phillips CM, Nerenberg PS, Drennan CL, Stultz CM. Physical Basis of Metal-Binding Specificity in Escherichia coli NikR. J Am Chem Soc. 2009;131(29):10220–10228. [PMC free article] [PubMed]
31. Onufriev A, Bashford D, Case DA. Exploring protein native states and large-scale conformational changes with a modified generalized born model. Proteins: Struct, Funct, Bioinf. 2004;55(2):383–394. [PubMed]
32. Zagrovic B, Pande V. Solvent viscosity dependence of the folding rate of a small protein: Distributed computing study. J Comput Chem. 2003;24(12):1432–1436. [PubMed]
33. Hornak V, Okur A, Rizzo RC, Simmerling C. HIV-1 protease flaps spontaneously close to the correct structure in simulations following manual placement of an inhibitor into the open state. J Am Chem Soc. 2006;128(9):2812–2813. [PMC free article] [PubMed]
34. Hornak V, Okur A, Rizzo RC, Simmerling C. HIV-1 protease flaps spontaneously open and reclose in molecular dynamics simulations. Proc Natl Acad Sci U S A. 2006;103(4):915–920. [PubMed]
35. Case DA, Cheatham TE, Darden T, Gohlke H, Luo R, Merz KM, Onufriev A, Simmerling C, Wang B, Woods RJ. The Amber biomolecular simulation programs. J Comput Chem. 2005;26(16):1668–1688. [PMC free article] [PubMed]
36. Hornak V, Abel R, Okur A, Strockbine B, Roitberg A, Simmerling C. Comparison of multiple amber force fields and development of improved protein backbone parameters. Proteins: Struct, Funct, Bioinf. 2006;65(3):712–725. [PubMed]
37. Sindhikara DJ, Kim S, Voter AF, Roitberg AE. Bad Seeds Sprout Perilous Dynamics: Stochastic Thermostat Induced Trajectory Synchronization in Biomolecules. J Chem Theory Comput. 2009;5(6):1624–1631. [PubMed]
38. Stote RH, Karplus M. Zinc-Binding in Proteins and Solution - a Simple but Accurate Nonbonded Representation. Proteins-Structure Function and Genetics. 1995;23(1):12–31. [PubMed]
39. Sakharov DV, Lim C. Zn protein simulations including charge transfer and local polarization effects. J Am Chem Soc. 2005;127(13):4921–4929. [PubMed]
40. Aqvist J, Warshel a. Computer-Simulation of the Initial Proton-Transfer Step in Human Carbonic Anhydrase-I. J Mol Biol. 1992;224(1):7–14. [PubMed]
41. Pang YP. Successful molecular dynamics simulation of two zinc complexes bridged by a hydroxide in phosphotriesterase using the cationic dummy atom method. Proteins-Structure Function and Genetics. 2001;45(3):183–189. [PubMed]
42. Pang YP, Xu K, El Yazal J, Prendergast FG. Successful molecular dynamics simulation of the zinc-bound farnesyltransferase using the cationic dummy atom approach. Protein Sci. 2000;9(10):1857–1865. [PubMed]
44. Frisch MJ, T GW, Schlegel HB, Scuseria GE, Robb MA, Cheeseman JR, Montgomery JA, Jr., Vreven T, Kudin KN, Burant JC, Millam JM, Iyengar SS, Tomasi J, Barone V, Mennucci B, Cossi M, Scalmani G, Rega N, Petersson GA, Nakatsuji H, Hada M, Ehara M, Toyota K, Fukuda R, Hasegawa J, Ishida M, Nakajima T, Honda Y, Kitao O, Nakai H, Klene M, Li X, Knox JE, Hratchian HP, Cross JB, Bakken V, Adamo C, Jaramillo J, Gomperts R, Stratmann RE, Yazyev O, Austin AJ, Cammi R, Pomelli C, Ochterski W, Ayala PY, Morokuma K, Voth GA, Salvador P, Dannenberg JJ, Zakrzewski VG, Dapprich S, Daniels AD, Strain MC, Farkas O, Malick DK, Rabuck AD, Raghavachari K, Foresman JB, Ortiz JV, Cui Q, Baboul AG, Clifford S, Cioslowski J, Stefanov BB, Liu G, Liashenko A, Piskorz P, Komaromi I, Martin RL, Fox DJ, Keith T, Al-Laham MA, Peng CY, Nanayakkara A, Challacombe M, Gill PMW, Johnson B, Chen W, Wong MW, Gonzalez C, Pople JA. Gaussian 03. Wallingford CT: 2004. C.02.
45. Riley KE, Merz KM. Assessment of density functional theory methods for the computation of heats of formation and ionization potentials of systems containing third row transition metals. J Phys Chem A. 2007;111(27):6044–6053. [PubMed]
46. Kannan KK, Ramanadham M, Jones TA. Structure, Refinement, and Function of Carbonic-Anhydrase Isozymes - Refinement of Human Carbonic Anhydrase-I. Ann N Y Acad Sci. 1984;429(Jun):49–60. [PubMed]
47. Besler BH, Merz KMJ, Kollman PA. Atomic Charges Derived from Semiempirical Methods. J Comput Chem. 1990;11(4):431–439.
48. Woods RJ, Khalil M, Pell W, Moffat SH, Smith JVH. Derivation of Net Atomic Charges from Molecular Electrostatic Potentials. J. Comp. Chem. 1990;11(3):297–310.
49. Chirlian LE, Francl MM. Atomic Charges Derived from Electrostatic Potentials: A Detailed Study. J. Comp. Chem. 1987;8(6):894–905.
50. Singh UC, Kollman PA. An approach to computing electrostatic charges for molecules. J. Comp. Chem. 1984;5(2):129–145.
51. Williams DE, Yan JM. Point-Charge Models for Molecules Derived from Least-Squares Fitting of the Electric Potential. Adv. Atomic Mol. Phys. 1988;23:87–130.
52. Cox SR, Williams DE. Representation of the Molecular Electrostatic Potential by a Net Atomic Charge Model. J Comput Chem. 1981;2(3):304–323.
53. Bayly CA, Cieplak P, Cornell WD, Kollman PA. A Well Behaved Elestrostatic Potential Based Method Using Charge Restraints for Deriving Atomic Charges: The RESP Model. J Phys Chem. 1993;97:10269–10280.
54. Floris FM, Persico M, Tani A, Tomasi J. Hydration Shell Structure of the Calcium-Ion from Simulations with Ab-Initio Effective Pair Potentials. Chem Phys Lett. 1994;227(1–2):126–132.
55. Bult CJ, White O, Olsen GJ, Zhou LX, Fleischmann RD, Sutton GG, Blake JA, FitzGerald LM, Clayton RA, Gocayne JD, Kerlavage AR, Dougherty BA, Tomb JF, Adams MD, Reich CI, Overbeek R, Kirkness EF, Weinstock KG, Merrick JM, Glodek A, Scott JL, Geoghagen NSM, Weidman JF, Fuhrmann JL, Nguyen D, Utterback TR, Kelley JM, Peterson JD, Sadow PW, Hanna MC, Cotton MD, Roberts KM, Hurst MA, Kaine BP, Borodovsky M, Klenk HP, Fraser CM, Smith HO, Woese CR, Venter JC. Complete genome sequence of the methanogenic archaeon, Methanococcus jannaschii. Science. 1996;273(5278):1058–1073. [PubMed]
56. Klenk HP, Clayton RA, Tomb JF, White O, Nelson KE, Ketchum KA, Dodson RJ, Gwinn M, Hickey EK, Peterson JD, Richardson DL, Kerlavage AR, Graham DE, Kyrpides NC, Fleischmann RD, Quackenbush J, Lee NH, Sutton GG, Gill S, Kirkness EF, Dougherty BA, McKenney K, Adams MD, Loftus B, Peterson S, Reich CI, McNeil LK, Badger JH, Glodek A, Zhou LX, Overbeek R, Gocayne JD, Weidman JF, McDonald L, Utterback T, Cotton MD, Spriggs T, Artiach P, Kaine BP, Sykes SM, Sadow PW, DAndrea KP, Bowman C, Fujii C, Garland SA, Mason TM, Olsen GJ, Fraser CM, Smith HO, Woese CR, Venter JC. The complete genome sequence of the hyperthermophilic, sulphate-reducing archaeon Archaeoglobus fulgidus. Nature. 1997;390(6658):364–370. [PubMed]
57. Tomb JF, White O, Kerlavage AR, Clayton RA, Sutton GG, Fleischmann RD, Ketchum KA, Klenk HP, Gill S, Dougherty BA, Nelson K, Quackenbush J, Zhou LX, Kirkness EF, Peterson S, Loftus B, Richardson D, Dodson R, Khalak HG, Glodek A, McKenney K, Fitzegerald LM, Lee N, Adams MD, Hickey EK, Berg DE, Gocayne JD, Utterback TR, Peterson JD, Kelley JM, Cotton MD, Weidman JM, Fujii C, Bowman C, Watthey L, Wallin E, Hayes WS, Borodovsky M, Karp PD, Smith HO, Fraser CM, Venter JC. The complete genome sequence of the gastric pathogen Helicobacter pylori (vol 388, pg 539, 1997) Nature. 1997;389(6649):412–412. [PubMed]
58. McClelland M, Sanderson KE, Spieth J, Clifton SW, Latreille P, Courtney L, Porwollik S, Ali J, Dante M, Du FY, Hou SF, Layman D, Leonard S, Nguyen C, Scott K, Holmes A, Grewal N, Mulvaney E, Ryan E, Sun H, Florea L, Miller W, Stoneking T, Nhan M, Waterston R, Wilson RK. Complete genome sequence of Salmonella enterica serovar typhimurium LT2. Nature. 2001;413(6858):852–856. [PubMed]
59. Deppenmeier U, Johann A, Hartsch T, Merkl R, Schmitz RA, Martinez-Arias R, Henne A, Wiezer A, Baumer S, Jacobi C, Bruggemann H, Lienard T, Christmann A, Bomeke M, Steckel S, Bhattacharyya A, Lykidis A, Overbeek R, Klenk HP, Gunsalus RP, Fritz HJ, Gottschalk G. The genome of Methanosarcina mazei: Evidence for lateral gene transfer between bacteria and archaea. J Mol Microbiol Biotechnol. 2002;4(4):453–461. [PubMed]
60. Eisen JA, Nelson KE, Paulsen IT, Heidelberg JF, Wu M, Dodson RJ, Deboy R, Gwinn ML, Nelson WC, Haft DH, Hickey EK, Peterson JD, Durkin AS, Kolonay JL, Yang F, Holt I, Umayam LA, Mason T, Brenner M, Shea TP, Parksey D, Nierman WC, Feldblyum TV, Hansen CL, Craven MB, Radune D, Vamathevan J, Khouri H, White O, Gruber TM, Ketchum KA, Venter JC, Tettelin H, Bryant DA, Fraser CM. The complete genome sequence of Chlorobium tepidum TLS, a photosynthetic, anaerobic, green-sulfur bacterium. Proc Natl Acad Sci U S A. 2002;99(14):9509–9514. [PubMed]
61. Galagan JE, Nusbaum C, Roy A, Endrizzi MG, Macdonald P, FitzHugh W, Calvo S, Engels R, Smirnov S, Atnoor D, Brown A, Allen N, Naylor J, Stange-Thomann N, DeArellano K, Johnson R, Linton L, McEwan P, McKernan K, Talamas J, Tirrell A, Ye WJ, Zimmer A, Barber RD, Cann I, Graham DE, Grahame DA, Guss AM, Hedderich R, Ingram-Smith C, Kuettner HC, Krzycki JA, Leigh JA, Li WX, Liu JF, Mukhopadhyay B, Reeve JN, Smith K, Springer TA, Umayam LA, White O, White RH, de Macario EC, Ferry JG, Jarrell KF, Jing H, Macario AJL, Paulsen I, Pritchett M, Sowers KR, Swanson RV, Zinder SH, Lander E, Metcalf WW, Birren B. The genome of M-acetivorans reveals extensive metabolic and physiological diversity. Genome Res. 2002;12(4):532–542. [PubMed]
62. Kaneko T, Nakamura Y, Sato S, Minamisawa K, Uchiumi T, Sasamoto S, Watanabe A, Idesawa K, Iriguchi M, Kawashima K, Kohara M, Matsumoto M, Shimpo S, Tsuruoka H, Wada T, Yamada M, Tabata S. Complete genomic sequence of nitrogen-fixing symbiotic bacterium Bradyrhizobium japonicum USDA 110. DNA Res. 2002;9(6):189–197. [PubMed]
63. Kawarabayasi Y, Sawada M, Horikawa H, Hino Y, Yamamoto S, Sekine M, Baba S, Kosugi H, Hosoyama A, Nagai Y, Sakai M, Ogura K, Otsuka R, Nakazawa H, Takamiya M, Ohfuku Y, Funahashi T, Tanaka T, Kudoh Y, Yamazaki J, Kushida N, Oguchi A, Aoki K, Kikuchi H. Complete sequence and gene organization of the genome of a hyper-thermophilic archaebacterium, Pyrococcus horikoshii OT3. DNA Res. 2002;5(2):55–76. [PubMed]
64. Paulsen IT, Seshadri R, Nelson KE, Eisen JA, Heidelberg JF, Read TD, Dodson RJ, Umayam L, Brinkac LM, Beanan MJ, Daugherty SC, Deboy RT, Durkin AS, Kolonay JF, Madupu R, Nelson WC, Ayodeji B, Kraul M, Shetty J, Malek J, Van Aken SE, Riedmuller S, Tettelin H, Gill SR, White O, Salzberg SL, Hoover DL, Lindler LE, Halling SM, Boyle SM, Fraser CM. The Brucella suis genome reveals fundamental similarities between animal and plant pathogens and symbionts. Proc Natl Acad Sci U S A. 2002;99(20):13148–13153. [PubMed]
65. Cohen GN, Barbe V, Flament D, Galperin M, Heilig R, Lecompte O, Poch O, Prieur D, Querellou J, Ripp R, Thierry JC, Van der Oost J, Weissenbach J, Zivanovic Y, Forterre P. An integrated analysis of the genome of the hyperthermophilic archaeon Pyrococcus abyssi. Mol Microbiol. 2003;47(6):1495–1512. [PubMed]
66. Mccammon JA, Gelin BR, Karplus M, Wolynes PG. Hinge-Bending Mode in Lysozyme. Nature. 1976;262(5566):325–326. [PubMed]
67. Mao B, Pear MR, Mccammon JA, Quiocho FA. Hinge-Bending in L-Arabinose-Binding Protein - the Venus-Flytrap Model. J Biol Chem. 1982;257(3):1131–1133. [PubMed]
68. Mao B, Mccammon JA. Structural Study of Hinge Bending in L-Arabinose-Binding Protein. J Biol Chem. 1984;259(8):4964–4970. [PubMed]