Computational protocol: Determinants of Cofactor Specificity for the Glucose-6-Phosphate Dehydrogenase from Escherichia coli: Simulation, Kinetics and Evolutionary Studies

Similar protocols

Protocol publication

[…] Homology models for the structure of EcG6PDH were created using the crystallographic structure of LmG6PDH in complex with NADP+ as the reference (PDB ID 1H9A). The alignment between the LmG6PDH and EcG6PDH sequences was obtained from the multiple sequences alignment generated as described in the section Phylogenetic Analysis of Bacterial G6PDHs. Ten three-dimensional models of EcG6PDH complexed to NADP+ were generated using Modeller 9.10 []. The consensus secondary structure prediction by Jpred 3 [] was used as a restraint during modeling. Subsequently, the models were ranked according to scoring functions implemented in Verify3D [] and Prosa2003 []. The best model was used for loop refinement (https://salilab.org/modeller/manual/node35.html) to generate an additional ten models that were also analyzed and scored. From the best final model, the structure of the EcG6PDH-NAD+ complex was obtained by removing the 2’-phosphate from the NADP+ complex. In order to take into account possible reorganization of the interactions in the binding site, 30,000 steps of conjugate gradient energy minimization were performed on the models. [...] Molecular Dynamics simulations of the model complexes of EcG6PDH with NADP+ or NAD+ were performed using NAMD (NAnoscale Molecular Dynamics) 2.8 [] with the Amber99SB-ILDN force field [] (which includes improved side-chain torsion potentials for the amino acids Isoleucine, Leucine, Aspartate and Asparagine) under explicit solvent, neutralized with Na+ ions in a box of TIP3P waters with a pad of 15 Å in all directions. The simulation was run using a time step of 1 fs, with a 12 Å cutoff for nonbonding interactions and a switching function from 10 Å for the Van der Waals and electrostatic terms. For the long-range interactions, the Particle Mesh Ewald method was applied. The simulation was carried out in the isothermal-isobaric (NPT) ensemble mode. The temperature was maintained constant through Langevin dynamics at 300 K, and the system was equilibrated until obtaining a constant RMSD for the backbone. Subsequently, the evolution of each system was simulated for 20 ns.In the analysis stage, the hydrogen bonds were quantified by considering an angle of greater than 120° between donor, hydrogen and acceptor atoms, and a distance below 3.5 Å between donor and acceptor. The binding free energies per residue were calculated by using the MM/PBSA (Molecular Mechanics Poisson—Boltzmann Surface Area) method implemented in AMBER 11 []. Snapshots for the calculations were taken every 50 ps, resulting in a total of 400 snapshots per trajectory. The ionic strength in molarity units was set to 0.1. Other options were set as default. In order to visualize the free energy at the surface of the protein, the values obtained for the backbone and side chains in the previous analysis were used to replace the β-factor column of each residue in the PDB file of the model, using Tcl/Tk scripts executed in VMD 1.9 []. [...] G6PDH K18A, R50A, K18A-R50A and K18T mutants were purified as described in []. The His-Tag was cut by TEV protease using a molar ratio of 1:50 (TEV:EcG6PDH) at 25°C, overnight. Throughout the different stages of purification, the enzyme activity was followed by the increment in absorbance at 340 nm, using the following solution: Tris-HCl 10mM pH 8.2, 1mM MgCl2, NADP+ 1mM and G6P 2mM. Protein concentrations were determined using the Bradford reagent (BioRad), with bovine serum albumin as the standard.Kinetic constants were determined at 25°C essentially as described by Olavarria et al 2012 []. Briefly, NAD(P)H formation was followed at 340 nm using a Synergy 2 spectrophotometer (BioTek) and 96-well plates, with a working volume of 300 μL. The substrates were added to the assay buffer immediately before loading the plate; reactions were triggered by addition of G6P. The assay was prepared using different concentrations of NADP+ or NAD+, at near-saturating concentrations of the co-substrate G6P (40 mM in the case of K18A and K18A/R50A, 50 mM for R50A; 15 mM for K18T; see ). NADP+, NAD+ and G6P were purchased from Sigma-Aldrich (96.5% purity for NAD+ [N7004], 95% for NADP+ [N5755] and premium quality for G6P [G7879] according to the vendor) solutions were pH neutralized, and concentrations were enzymatically titred.To obtain a good approximation of the initial rate, slopes were taken before 5% of initial substrate was consumed. Three independent measurements were taken and the standard deviation was calculated. Data was adjusted to the Michaelis-Menten equation to obtain KM and kcat, by using SigmaPlot (Systat Software, San Jose, CA). [...] In order to generate a phylogenetic tree of bacterial G6PDHs, we first searched the Swiss-Prot and TrEMBL (UniProtKB) databases for sequences of enzymes that have been kinetically characterized using NAD+ and NADP+ as cofactor. For 6 cases, we were able to find 1 G6PDH sequence among the bacterial species (or strains) for which kinetic characterization had been reported. For the remainder, either a specific strain was not found in UniprotKB (in other words, the sequence was annotated only up to the level of species), or several genes could be associated to the bacterial species or strain from which the G6PDH was characterized (authors did not report the specific isoform that was studied). We retrieved a total of 53 G6PDH sequences. After alignment of these sequences using ClustalX [] we calculated a Neighbor Joining tree and focused on clusters of sequences from single organisms. Within each cluster, several sequences showed the same residue at positions 18 and 50 (EcG6PDH numbering), and so only one representative was selected. After this step, 17 sequences remained, including representatives from the phyla Proteobacteria (α, β & γ), Firmicutes, Thermotogae and Aquificae. Four additional G6PDH sequences were included: those from Mycobacterium avium (Actinobacteria), whose structure is known (PDB ID 4LGV), from Borreliaburgdorferi (Spirochaetes), Synechocisits (Cyanobacteria) and Chlamydophila pneumoniae (Chlamydiae). The latter were included in order to enrich the number of analyzed phyla. In a preliminary tree we observed that G6PDH from M. avium coalesced with three proteobacterial sequences (Gluconacetobacter hansenii b and c, and Pseudomonas fluorescens b) to form a separate cluster. To better describe sequence conservation within this group, we searched the databases for sequences with less than 70% identity (to avoid redundancy). We ended up working with 31 bacterial G6PDH sequences (). The structural alignment of the G6PDHs from Leuconostoc mesenteroides (L. mesenteroides) and Mycobacterium avium (M. avium) was used as an alignment profile in ClustalX to which the remainder of the sequences were subsequently aligned. This alignment was used to construct a Bayesian tree, using MrBayes [] in a calculation of 5x106 generations. Trees were sampled every 500 generations and the analysis was performed discarding the first 2500 trees. […]

Pipeline specifications

Software tools MODELLER, JPred, VERIFY3D, NAMD, AMBER, VMD, SigmaPlot, Clustal W, MrBayes
Databases UniProt UniProtKB
Applications Miscellaneous, Phylogenetics, Protein structure analysis
Organisms Escherichia coli
Chemicals Alanine, NAD, NADP, Glucose-6-Phosphate