Study of the docking of competitive inhibitors at a model of tyrosinase active site: insights from joint broken-symmetry/Spin-Flip DFT computations and ELF topological analysis

Following our previous study (Piquemal et al., New J. Chem., 2003, 27, 909), we present here a DFT study of the inhibition of the Tyrosinase enzyme. Broken-symmetry DFT computations are supplemented with Spin-Flip TD-DFT calculations, which, for the first time, are applied to such a dicopper enzyme. The chosen biomimetic model encompasses a dioxygen molecule, two Cu(II) cations, and six imidazole rings. The docking energy of a natural substrate, namely phenolate, together with those of several inhibitor and non-inhibitor compounds, are reported and show the ability of the model to rank the most potent inhibitors in agreement with experimental data. With respect to broken-symmetry calculations, the Spin-Flip TD-DFT approach reinforces the possibility for theory to point out potent inhibitors: the need for the deprotonation of the substrates, natural or inhibitors, is now clearly established. Moreover, Electron Localization Function (ELF) topological analysis computations are used to deeply track the particular electronic distribution of the Cu-O-Cu three-center bonds involved in the enzymatic Cu(2)O(2) metallic core (Piquemal and Pilmé, J. Mol. Struct.: Theochem, 2006, 77, 764). It is shown that such bonds exhibit very resilient out-of-plane density expansions that play a key role in docking interactions: their 3D-orientation could be the topological electronic signature of oxygen activation within such systems.

Simple Formulas for Improved Point-Charge Electrostatics in Classical Force Fields and Hybrid Quantum Mechanical/Molecular Mechanical Embedding

We present a simple damping scheme for point-charge electrostatics that could be used directly in classical force fields. The approach acts at the charge (or monopole) level only and allows the inclusion of short-range electrostatic penetration effects at a very low cost. Results are compared with density functional theory Coulomb intermolecular interaction energies and with several other methods such as distributed multipoles, damped distributed multipoles, and transferable Hermite-Gaussian densities. Realistic trends in the interactions are observed for atom-centered Mertz-Kollman corrected point-charge distributions. The approach allows increasing the selectivity of parameters in the case of metal complexes. In addition, two QM/MM calculations are presented where the damping function is employed to include the MM atoms located at the QM/MM boundary. The first calculation corresponds to the gas-phase proton transfer of aspartic acid through water and the second is the first step of the reaction catalyzed by the 4-oxalocrotonate tautomerase (4OT) enzyme. First, improved agreement is observed when using the damping approach compared with the conventional excluded charge method or when including all charges in the calculation. Second, in the case of 4OT, the damped charge approach is in agreement with previous calculations, whereas including all charges gives a significantly higher energy barrier. In both cases, no reparameterization of the van der Waals part of the force field was performed.

An all-atom solution-equilibrated model for human extrinsic blood coagulation complex (sTF-VIIa-Xa): a protein-protein docking and molecular dynamics refinement study

Tissue factor (TF)-bound factor (F)VIIa plays a critical role in activating FX, an event that rapidly results in blood coagulation. Despite recent advances in the structural information about soluble TF (sTF)-bound VIIa and Xa individually, the atomic details of the ternary complex are not known. As part of our long-term goal to provide a structural understanding of the extrinsic blood coagulation pathway, we built an all atom solution-equilibrated model of the human sTF-VIIa-Xa ternary complex using protein-protein docking and molecular dynamics (MD) simulations. The starting structural coordinates of sTF-VIIa and Xa were derived from dynamically equilibrated solution structures. Due to the flexible nature of the light-chain of the Xa molecule, a three-stage docking approach was employed in which SP (Arg195-Lys448)/EGF2 (Arg86-Arg139), EGF1 (Asp46-Thr85) and GLA (Ala1-Lys45) domains were docked in a sequential manner. The rigid-body docking approach of the FTDOCK method in conjunction with filtering based on biochemical knowledge from experimental site-specific mutagenesis studies provided the strategy. The best complex obtained from the docking experiments was further refined using MD simulations for 3 ns in explicit water. In addition to explaining most of the known experimental site-specific mutagenesis data pertaining to sTF-VIIa, our model also characterizes likely enzyme-binding exosites on FVIIa and Xa that may be involved in the ternary complex formation. According to the equilibrated model, the 140s loop of VIIa serves as the key recognition motif for complex formation. Stable interactions occur between the FVIIa 140s loop and the FXa -strand B2 region near the sodium-binding domain, the 160 s loop and the N-terminal activation loop regions. The helical-hydrophobic stack region that connects the GLA and EGF1 domains of VIIa and Xa appears to play a potential role in the membrane binding region of the ternary complex. The proposed model may serve as a reasonable structural basis for understanding the exosite-mediated substrate recognition of sTF-VIIa and to advance understanding of the TFPI-mediated regulatory pathway of the extrinsic blood coagulation cascade.

Gene selection for sample classification based on gene expression data: study of sensitivity to choice of parameters of the GA/KNN method

MOTIVATION

We recently introduced a multivariate approach that selects a subset of predictive genes jointly for sample classification based on expression data. We tested the algorithm on colon and leukemia data sets. As an extension to our earlier work, we systematically examine the sensitivity, reproducibility and stability of gene selection/sample classification to the choice of parameters of the algorithm.

METHODS

Our approach combines a Genetic Algorithm (GA) and the k-Nearest Neighbor (KNN) method to identify genes that can jointly discriminate between different classes of samples (e.g. normal versus tumor). The GA/KNN method is a stochastic supervised pattern recognition method. The genes identified are subsequently used to classify independent test set samples.

RESULTS

The GA/KNN method is capable of selecting a subset of predictive genes from a large noisy data set for sample classification. It is a multivariate approach that can capture the correlated structure in the data. We find that for a given data set gene selection is highly repeatable in independent runs using the GA/KNN method. In general, however, gene selection may be less robust than classification.

AVAILABILITY

The method is available at http://dir.niehs.nih.gov/microarray/datamining

CONTACT

LI3@niehs.nih.gov

Modeling human zymogen factor IX

Modern theoretical techniques are employed to provide complete three dimensional structure for the zymogen and activated forms of human coagulation factors IX and IXa. These structures are fully calcium bound and equilibrated in an electrically neutral aqueous environment. The relationship of structure to mutational data is examined. We find that a substantial relative orientational change of the catalytic domain occurs on activation. Also, we find that the electrostatistically dipolar nature of the catalytic domain is substantially modified upon activation, with cleavage of the negatively charged activation peptide leaving behind a largely hydrophobic face in factor IXa. While the backbone atoms of the catalytic residues have little relative movement, nearby loops are found that do move. The presence or absence of these changes likely defines specificity.

Biased distribution of inverted and direct Alus in the human genome: implications for insertion, exclusion, and genome stability

Alu sequences, the most abundant class of large dispersed DNA repeats in human chromosomes, contribute to human genome dynamics. Recently we reported that long inverted repeats, including human Alus, can be strong initiators of genetic change in yeast. We proposed that the potential for interactions between adjacent, closely related Alus would influence their stability and this would be reflected in their distribution. We have undertaken an extensive computational analysis of all Alus (the database is at http://dir.niehs.nih.gov/ALU) to better understand their distribution and circumstances under which Alu sequences might affect genome stability. Alus separated by <650 bp were categorized according to orientation, length of regions sharing high sequence identity, distance between highly identical regions, and extent of sequence identity. Nearly 50% of all Alu pairs have long alignable regions (>275 bp), corresponding to nearly full-length Alus, regardless of orientation. There are dramatic differences in the distributions and character of Alu pairs with closely spaced, nearly identical regions. For Alu pairs that are directly repetitive, approximately 30% have highly identical regions separated by <20 bp, but only when the alignments correspond to near full-size or half-size Alus. The opposite is found for the distribution of inverted repeats: Alu pairs with aligned regions separated by <20 bp are rare. Furthermore, closely spaced direct and inverted Alus differ in their truncation patterns, suggesting differences in the mechanisms of insertion. At larger distances, the direct and inverted Alu pairs have similar distributions. We propose that sequence identity, orientation, and distance are important factors determining insertion of adjacent Alus, the frequency and spectrum of Alu-associated changes in the genome, and the contribution of Alu pairs to genome instability. Based on results in model systems and the present analysis, closely spaced inverted Alu pairs with long regions of alignment are likely at-risk motifs (ARMs) for genome instability.

Heparan sulfate biosynthesis: a theoretical study of the initial sulfation step by N-deacetylase/N-sulfotransferase

Heparan sulfate N-deacetylase/N-sulfotransferase (NDST) catalyzes the deacetylation and sulfation of N-acetyl-D-glucosamine residues of heparan sulfate, a key step in its biosynthesis. Recent crystallographic and mutational studies have identified several potentially catalytic residues of the sulfotransferase domain of this enzyme (, J. Biol. Chem. 274:10673-10676). We have used the x-ray crystal structure of heparan sulfate N-sulfotransferase with 3'-phosphoadenosine 5'-phosphate to build a solution model with cofactor 3'-phosphoadenosine 5'-phosphosulfate (PAPS) and a model heparan sulfate ligand bound, and subsequently performed a 2-ns dynamics solution simulation. The simulation results confirm the importance of residues Glu(642), Lys(614), and Lys(833), with the possible involvement of Thr(617) and Thr(618), in binding PAPS. Additionally, Lys(676) is found in close proximity to the reaction site in our solvated structure. This study illustrates for the first time the possible involvement of water in the catalysis. Three water molecules were found in the binding site, where they are coordinated to PAPS, heparan sulfate, and the catalytic residues.

Modeling zymogen protein C

A solution structure for the complete zymogen form of human coagulation protein C is modeled. The initial core structure is based on the x-ray crystallographic structure of the gamma-carboxyglutamic acid (Gla)-domainless activated form. The Gla domain (residues 1-48) is modeled from the x-ray crystal coordinates of the factor VII(a)/tissue factor complex and oriented with the epidermal growth factor-1 domain to yield an initial orientation consistent with the x-ray crystal structure of porcine factor IX(a). The missing C-terminal residues in the light chain (residues 147-157) and the activation peptide residues 158-169 were introduced using homology modeling so that the activation peptide residues directly interact with the residues in the calcium binding loop. Molecular dynamics simulations (Amber-particle-mesh-Ewald) are used to obtain the complete calcium-complexed solution structure. The individual domain structures of protein C in solution are largely unaffected by solvation, whereas the Gla-epidermal growth factor-1 orientation evolves to a form different from both factors VII(a) and IX(a). The solution structure of the zymogen protein C is compared with the crystal structures of the existing zymogen serine proteases: chymotrypsinogen, proproteinase, and prethrombin-2. Calculated electrostatic potential surfaces support the involvement of the serine protease calcium ion binding loop in providing a suitable electrostatic environment around the scissile bond for II(a)/thrombomodulin interaction.

Heparan/chondroitin sulfate biosynthesis. Structure and mechanism of human glucuronyltransferase I

Human beta1,3-glucuronyltransferase I (GlcAT-I) is a central enzyme in the initial steps of proteoglycan synthesis. GlcAT-I transfers a glucuronic acid moiety from the uridine diphosphate-glucuronic acid (UDP-GlcUA) to the common linkage region trisaccharide Gal beta 1-3Gal beta 1-4Xyl covalently bound to a Ser residue at the glycosaminylglycan attachment site of proteoglycans. We have now determined the crystal structure of GlcAT-1 at 2.3 A in the presence of the donor substrate product UDP, the catalytic Mn(2+) ion, and the acceptor substrate analog Gal beta 1-3Gal beta 1-4Xyl. The enzyme is a alpha/beta protein with two subdomains that constitute the donor and acceptor substrate binding site. The active site residues lie in a cleft extending across both subdomains in which the trisaccharide molecule is oriented perpendicular to the UDP. Residues Glu(227), Asp(252), and Glu(281) dictate the binding orientation of the terminal Gal-2 moiety. Residue Glu(281) is in position to function as a catalytic base by deprotonating the incoming 3-hydroxyl group of the acceptor. The conserved DXD motif (Asp(194), Asp(195), Asp(196)) has direct interaction with the ribose of the UDP molecule as well as with the Mn(2+) ion. The key residues involved in substrate binding and catalysis are conserved in the glucuronyltransferase family as well as other glycosyltransferases.

Molecular dynamics simulations of the d(CCAACGTTGG)(2) decamer: influence of the crystal environment

Molecular dynamics (MD) simulations of the DNA duplex d(CCAACGTTGG)(2) were used to study the relationship between DNA sequence and structure. Two crystal simulations were carried out; one consisted of one unit cell containing two duplexes, and the other of two unit cells containing four duplexes. Two solution simulations were also carried out, one starting from canonical B-DNA and the other starting from the crystal structure. For many helicoidal parameters, the results from the crystal and solution simulations were essentially identical. However, for other parameters, in particular, alpha, gamma, delta, (epsilon - zeta), phase, and helical twist, differences between crystal and solution simulations were apparent. Notably, during crystal simulations, values of helical twist remained comparable to those in the crystal structure, to include the sequence-dependent differences among base steps, in which values ranged from 20 degrees to 50 degrees per base step. However, in the solution simulations, not only did the average values of helical twist decrease to approximately 30 degrees per base step, but every base step was approximately 30 degrees, suggesting that the sequence-dependent information may be lost. This study reveals that MD simulations of the crystal environment complement solution simulations in validating the applicability of MD to the analysis of DNA structure.

Molecular dynamics simulations of biomolecules: long-range electrostatic effects

Current computer simulations of biomolecules typically make use of classical molecular dynamics methods, as a very large number (tens to hundreds of thousands) of atoms are involved over timescales of many nanoseconds. The methodology for treating short-range bonded and van der Waals interactions has matured. However, long-range electrostatic interactions still represent a bottleneck in simulations. In this article, we introduce the basic issues for an accurate representation of the relevant electrostatic interactions. In spite of the huge computational time demanded by most biomolecular systems, it is no longer necessary to resort to uncontrolled approximations such as the use of cutoffs. In particular, we discuss the Ewald summation methods, the fast particle mesh methods, and the fast multipole methods. We also review recent efforts to understand the role of boundary conditions in systems with long-range interactions, and conclude with a short perspective on future trends.

Residues in the alphaH and alphaI helices of the HIV-1 reverse transcriptase thumb subdomain required for the specificity of RNase H-catalyzed removal of the polypurine tract primer

During retrovirus replication, reverse transcriptase (RT) must specifically interact with the polypurine tract (PPT) to generate and subsequently remove the RNA primer for plus-strand DNA synthesis. We have investigated the role that human immunodeficiency virus-1 RT residues in the alphaH and alphaI helices in the thumb subdomain play in specific RNase H cleavage at the 3'-end of the PPT; an in vitro assay modeling the primer removal step was used. Analysis of alanine-scanning mutants revealed that a subgroup exhibits an unusual phenotype in which the PPT is cleaved up to seven bases from its 3'-end. Further analysis of alphaH mutants (G262A, K263A, N265A, and W266A) with changes in residues in or near a structural motif known as the minor groove binding track showed that the RNase H activity of these mutants is more dramatically affected with PPT substrates than with non-PPT substrates. Vertical scan mutants at position 266 were all defective in specific RNase H cleavage, consistent with conservation of tryptophan at this position among lentiviral RTs. Our results indicate that residues in the thumb subdomain and the minor groove binding track in particular, are crucial for unique interactions between RT and the PPT required for correct positioning and precise RNase H cleavage.

Probing the structural changes in the light chain of human coagulation factor VIIa due to tissue factor association

The crystallographic structure of human coagulation factor VIIa/tissue factor complex bound with calcium ions was used to model the solution structure of the light chain of factor VIIa (residues 1-142) in the absence of tissue factor. The Amber force field in conjunction with the particle mesh Ewald summation method to accommodate long-range electrostatic interactions was used in the trajectory calculations. The estimated TF-free solution structure was then compared with the crystal structure of factor VIIa/tissue factor complex to estimate the restructuring of factor VIIa due to tissue factor binding. The solution structure of the light chain of factor VIIa in the absence of tissue factor is predicted to be an extended domain structure similar to that of the tissue factor-bound crystal. Removal of the EGF1-bound calcium ion is shown by simulation to lead to minor structural changes within the EGF1 domain, but also leads to substantial relative reorientation of the Gla and EGF1 domains.

An atomic model for the pleated beta-sheet structure of Abeta amyloid protofilaments

Synchrotron x-ray studies on amyloid fibrils have suggested that the stacked pleated beta-sheets are twisted so that a repeating unit of 24 beta-strands forms a helical turn around the fibril axis (. J. Mol. Biol. 273:729-739). Based on this morphological study, we have constructed an atomic model for the twisted pleated beta-sheet of human Abeta amyloid protofilament. In the model, 48 monomers of Abeta 12-42 stack (four per layer) to form a helical turn of beta-sheet. Each monomer is in an antiparallel beta-sheet conformation with a turn located at residues 25-28. Residues 17-21 and 31-36 form a hydrophobic core along the fibril axis. The hydrophobic core should play a critical role in initializing Abeta aggregation and in stabilizing the aggregates. The model was tested using molecular dynamics simulations in explicit aqueous solution, with the particle mesh Ewald (PME) method employed to accommodate long-range electrostatic forces. Based on the molecular dynamics simulations, we hypothesize that an isolated protofilament, if it exists, may not be twisted, as it appears to be when in the fibril environment. The twisted nature of the protofilaments in amyloid fibrils is likely the result of stabilizing packing interactions of the protofilaments. The model also provides a binding mode for Congo red on Abeta amyloid fibrils. The model may be useful for the design of Abeta aggregation inhibitors.

Vertical-scanning mutagenesis of a critical tryptophan in the minor groove binding track of HIV-1 reverse transcriptase. Molecular nature of polymerase-nucleic acid interactions

While sequence-specific DNA-binding proteins interact predominantly in the DNA major groove, DNA polymerases bind DNA through interactions in the minor groove that are sequence nonspecific. Through functional analyses of alanine-substituted mutant enzymes that were guided by molecular dynamics modeling of the human immunodeficiency virus type 1-reverse transcriptase and DNA complex, we previously identified a structural element in reverse transcriptase, the minor groove binding track (MGBT). The MGBT is comprised of five residues (Ile94, Gln258, Gly262, Trp266, and Gln269) which interact 2-6 base pairs upstream from the polymerase active site in the DNA minor groove and are important in DNA binding, processivity, and frameshift fidelity. These residues do not contribute equally; functional analysis of alanine mutants suggests that Trp266 contributes the most to binding. To define the molecular interactions between Trp266 and the DNA minor groove, we have analyzed the properties of eight mutants, each with an alternate side chain at this position. A refined molecular dynamics model was used to calculate relative binding free energies based on apolar surface area buried upon complex formation. In general, there was a strong correlation between the relative calculated binding free energies for the alternate residue 266 side chains and the magnitude of the change in the properties which reflect template-primer interactions (template-primer dissociation rate constant, Ki,AZTTP, processivity, and frameshift fidelity). This correlation suggests that hydrophobic interactions make a major contribution to the stability of the polymerase-DNA complex. Additionally, tyrosine and arginine substitutions resulted in mutant enzymes with DNA binding properties better than predicted by buried surface area alone, suggesting that hydrogen bonding could also play a role in DNA binding at this position.

Trans-cis isomerization of proline 22 in bovine prothrombin fragment 1: a surprising result of structural characterization

The calcium ion-mediated interaction of bovine prothrombin (BF1) with negatively charged phospholipid membranes is assumed to be largely via the Gla domain of BF1 with the fold of the Gla domain essential for binding. It has been reported that Pro22 undergoes classical trans to cis isomerization in the presence of calcium ions with the cis conformation of Pro22 of BF1 responsible for membrane binding [Evans, T. C., Jr., and Nelsestuen, G. L. (1996) Biochemistry 35, 8210-8215]. However, Pro22 was found to be in the trans conformation in the crystal structure of BF1. In the present work, we have used molecular dynamics simulations to investigate the relative importance of the two conformations of Pro22 to the structural and dynamical properties of BF1. The initial trans conformation of Pro22 in BF1 was slowly converted to cis-Pro22 using constrained dynamics. The second-generation AMBER force field in conjunction with the particle mesh Ewald method to accommodate long-range interaction was employed in the trajectory calculations. Comparison of the BF1(trans-Pro22) and BF1(cis-Pro22) equilibrated structures reveals surprisingly that the overall structural changes associated with the trans-cis isomerization is minimal and only minor modifications to the hydrogen bond network and the network of N-terminus Ala1 take place. The calculated electrostatic potential energy surfaces of the two protein structures also appear to be very similar, indicating the near equality of the local interaction site environments in the protein prior to lipid binding.

Reciprocal size-effect relationship of the key residues in determining regio- and stereospecificities of DHEA hydroxylase activity in P450 2a5

Collectively, the P450 2a4/2a5 system hyrdoxylates DHEA in at least three positions (7alpha, 7beta, and 2alpha). An individual P450, however, exhibits high specificity to one of these products. Using site-directed mutagenesis of mP450 2a5 from the wild mouse Mus minutoides and bacterial expression, we have associated the function of residues 117, 209, and 481 with the respective specificity observed in each P450. Ala at position 117 determines the 7beta-hydroxylase activity, whereas Val at this position defines the 2alpha-hydroxylase activity. Leu at position 209 is essential for high DHEA 7alpha-hydroxylase activity. The substitutions of residue 481 with various hydrophobic amino acids elicited a profound alteration of the specific hydroxylation rates, but did not influence the regio- and stereospecificities at either of the three positions of DHEA. The alterations caused by residue 481 also depended on the residue identity at position 117 or 209. The results indicate that the sizes of several key residues obey a concerted reciprocal relationship whereby the substrate pocket of the P450s adjusts to accommodate DHEA. A limited molecular modeling study successfully correlates DHEA binding to experimental DHEA hydroxylase activities for a series of mutants at key positions.

Refinement of the NMR solution structure of the gamma-carboxyglutamic acid domain of coagulation factor IX using molecular dynamics simulation with initial Ca2+ positions determined by a genetic algorithm

A genetic algorithm (GA) successfully identified the calcium positions in the crystal structure of bovine prothrombin fragment 1 bound with calcium ions (bf1/Ca). The same protocol was then used to determine the calcium positions in a closely related fragment, the Gla domain of coagulation factor IX, the structure of which had previously been determined by NMR spectroscopy in the presence of calcium ions. The most frequently occurring low-energy structure found by GA was used as the starting structure for a molecular dynamics refinement. The molecular dynamics simulation was performed using explicit water and the Particle-Mesh Ewald method to accommodate the long-range electrostatic forces. While the overall conformation of the NMR structure was preserved, significant refinement is apparent when comparing the simulation average structure with its NMR precursor in terms of the N-terminal (Tyr1-N) network, the total number of hydrogen bonds, the calcium ion coordinations, and the compactness of the structure. It is likely that the placement of calcium ions in the protein is critical for refinement. The calcium ions apparently induce structural changes during the course of the simulation that result in a more compact structure.

Structural flexibility and functional versatility of mammalian P450 enzymes

P450 enzymes have evolved into a large superfamily that displays great diversity in substrate and product specificities by fixing the natural amino acid substitutions with high frequency. Site-directed mutagenesis has been used to correlate the substitutions with the diverse specificities in various P450s. As a result, the common residues that determine the specificities of various mammalian P450s have been identified and aligned to the corresponding residues in the substrate-heme pocket of the 3-dimensional structures of bacterial P450s. The substrate-heme pocket appears to be structurally variable so that only a minor substitution (Ala -> -> Val, for example) at the critical positions is enough to define the altered specificity. Thus, the structural variability of the P450s provides the inherent versatility in acquiring a novel activity. Recent mutational studies indicate that the side chain size is the major determining factor of specificity, outweighing other factors such as polarity. Further understanding of the paradoxical characteristics observed may provide us with the underlying principles that determine P450 activities, and may lead to the ability to predict P450 activities based on the types of key amino acid residues.

Selected new developments in computational chemistry

Molecular dynamics is a general technique for simulating the time-dependent properties of molecules and their environments. Quantum mechanics, as applied to molecules or clusters of molecules, provides a prescription for predicting properties exactly (in principle). It is reasonable to expect that both will have a profound effect on our understanding of environmental chemistry in the future. In this review, we consider several recent advances and applications in computational chemistry.

Structural flexibility and functional versatility of cytochrome P450 and rapid evolution

P450 represents a large group of heme-thiolate enzymes that exhibit remarkably diverse activities for the metabolism of numerous endogenous and exogenous chemicals. Recent site-directed mutagenesis studies indicate that a single mutation at any of the key residues can be enough to alter the substrate and/or product specificities in the P450 activities. Molecular modeling predicts that these key residues are located within the substrate heme pocket. Structural elements involved in diversifying P450 activity appear to correspond to the B' helix, the F helix and the F/G interhelical loop in the bacterial P450s. Structures represented by these regions are extremely variable despite the fact that the core of the P450 substrate pocket is well conserved. A mutation within these regions may result in a significant geometrical alteration of the pocket and lead to diversify the P450 activity. Phylogenetical analysis shows a relatively high rate of nonsynonymous substitution within these substrate binding regions. The functional versatility of P450 can thus be largely accounted for in terms of pocket change brought about by rapid mutations.

The roles of individual amino acids in altering substrate specificity of the P450 2a4/2a5 enzymes

A single amino acid substitution is sufficient to alter substrate specificity of P450 enzymes. Mouse P450 2a5, for example, has its substrate specificity converted from coumarin 7- to testosterone 15 alpha-hydroxylase activity by the substitution of Phe at position 209 to Leu. Furthermore, placing Asn at this position confers a novel corticosterone 15 alpha-hydroxylase activity to this P450. Recent site-directed mutational studies show the presence of the topologically common residues, each of which can determine the specificities of various mammalian P450s. For instance, residue 209 (in 2a5) corresponds to a residue at position 206 in rat P4502B1 that regulates its steroid hydroxylase activity. High substrate specificity often observed in an individual P450, therefore, can be determined and altered by the identities of a few critical residues. The structural flexibility of the substrate-heme pocket may also provide P450 enzymes with the ability to display a broad range of substrate specificities. Understanding the underlying principles whereby the flexible pocket determines P450 activities may lead us to the prediction of P450 activities based on the identities of key amino acid residues.

Reduced frameshift fidelity and processivity of HIV-1 reverse transcriptase mutants containing alanine substitutions in helix H of the thumb subdomain

We have analyzed two human immunodeficiency virus (HIV-1) reverse transcriptase mutants of helix H in the thumb subdomain suggested by x-ray crystallography to interact with the primer strand of the template-primer. These enzymes, G262A and W266A, were previously shown to have greatly elevated dissociation rate constants for template-primer and to be much less sensitive to inhibition by 3'-azidodeoxythymidine 5'-triphosphate. Here we describe their processivity and error specificity. The results reveal that: (i) both enzymes have reduced processivity and lower fidelity for template-primer slippage errors, (ii) they differ from each other in sequence-dependent termination of processive synthesis and in error specificity, and (iii) the magnitude of the mutator effect relative to wild-type enzyme for deletions in homopolymeric sequences decreases as the length of the run increases. Thus amino acid substitutions in a subdomain thought to interact with the duplex template-primer confer a strand slippage mutator phenotype to a replicative DNA polymerase. This suggests that interactions between specific amino acids and the primer stem at positions well removed from the active site are critical determinants of processivity and fidelity. These effects, obtained in aqueous solution during catalytic cycling, are consistent with and support the existing crystallographic structural model.

Altering the regiospecificity of androstenedione hydroxylase activity in P450s 2a-4/5 by a mutation of the residue at position 481

Mouse P450 2a-5 (coumarin 7-hydroxylase) acquires androstenedione (AD) hydroxylase activity by substituting Phe at position 209 with Asn. However, this mutant P450 2a-5 (F209N) and the corresponding mutant P450 2a-4 (L209N) exhibit different regiospecificites of androstenedione (AD) hydroxylase activity. While the former mutant catalyzes both AD 15 alpha- and 7 alpha-hydroxylase activities at similar rates, the latter mutant maintains the original high specificity of AD 15 alpha-hydroxylase activity. The AD hydroxylase activities in chimeric enzymes of the mutants L209N and F209N show that the regiospecificites are determined by the carboxy-terminal halves of the P450 molecules. Mutations at each of the four different residues within the carboxy-terminal halves indicate that the differences in regiospecificity are determined by the Val/Ala mutation at position 481. As the size of the hydrophobic amino acid at position 481 becomes larger (Ala < Val < Ile), the regiospecificities toward the C15 position of the AD molecule are dramatically increased. The regiospecificity is also increased by placing positively-charged Arg at position 481, although the remaining 15 alpha-hydroxylase activity in this mutant is considerably lower than the other P450s. The results indicate that the size of the residue at position 481 is a key factor in regulating the regiospecificity of AD hydroxylase activity in the P450s. Modeling AD in the substrate-heme pocket of bacterial P450 101A provided further support that residue 481 may reside near the steroid molecule so as to possibly affect the AD hydroxylase activity.

Multiple steroid-binding orientations: alteration of regiospecificity of dehydroepiandrosterone 2- and 7-hydroxylase activities of cytochrome P-450 2a-5 by mutation of residue 209

The mutation of Ala-117 to Val conferred dehydroepiandrosterone (DHEA) hydroxylase activity on cytochrome P-450 2a-4, with the production of both 2 alpha- and 7 alpha-hydroxyDHEA at similar rates. P-450 2a-5 which has Val at position 117, acquired high DHEA hydroxylase activity by mutation of Phe-209. Mutant F209L of P-450 2a-5 exhibited strong regiospecificity at the 2-position of the DHEA molecule with the production of 2 alpha-hydroxy DHEA as the major metabolite. On the other hand, mutant F209V of P-450 2a-5 showed the 7-position to be the major hydroxylation site, 7 beta-hydroxyDHEA and 7 alpha-OHDHEA being produced. Therefore the regiospecificity of DHEA hydroxylase activity of P-450 2a-5 is altered between the 2- and 7-position depending on the amino acid at position 209. Modelling of the DHEA molecule in the pocket of bacterial P-450cam showed that the steroid can be accommodated in at least two orientations for which the 2- or 7- position is near the sixth axial position of the haem. Moreover, these two orientations, which are of similar energy, can be interconverted by a 180 degrees rotation of the steroid molecule around its long axis. These results support the hypothesis that the steroid molecule in the pocket is in dynamic equilibrium with multiple binding orientations and that the equilibrium is apparently determined by a few critical residues including those at positions 117 and 209.

Atomic-level accuracy in simulations of large protein crystals

Proper treatment of long-range Coulombic forces presents a major obstacle to providing realistic molecular dynamics simulations of macromolecules. Traditional approximations made to lessen computational cost ultimately lead to unrealistic behavior. The particle mesh Ewald method accommodates long-range Coulombic forces accurately and efficiently by use of fast Fourier transform techniques. We report a 1-ns simulation of bovine pancreatic trypsin inhibitor in a crystal unit cell using the particle mesh Ewald methodology. We find an rms backbone deviation from the x-ray structure (0.33 A) that is lower than that observed between bovine pancreatic trypsin inhibitor in different crystal forms and much lower than those of previous simulations. These results bridge the gap between structures obtained from molecular simulation and those from experiment.

Inherent versatility of P-450 oxygenase. Conferring dehydroepiandrosterone hydroxylase activity to P-450 2a-4 by a single amino acid mutation at position 117

Mouse steroid 15 alpha-hydroxylase P-450 2a-4 is restricted in its substrate specificity to the delta 4, 3-ketone steroids such as androstenedione. As a result, the P-450 exhibits little hydroxylase activity toward delta 5, 3-hydroxysteroids including dehydroepiandrosterone (DHEA). A single amino acid mutation of Ala at position 117 to Val, however, is enough to confer a high DHEA hydroxylase activity to P-450 2a-4 with 7 alpha-OH DHEA as one of the two major hydroxylated metabolites. Mouse coumarin 7-hydroxylase P-450 2a-5 contains Val at position 117, but it exhibits very low DHEA hydroxylase activity. P-450 2a-5 acquires high DHEA hydroxylase activity, however, by a mutation of Phe-209 to Asn. Moreover, the mutant P-450 2a-5 loses its activity when Val is replaced by Ala at position 117. The residue at position 117, therefore, plays the principal role in the determination of the DHEA hydroxylase activity of the P-450s. Conversely, mutations at residue 117 have little effect on the androstenedione hydroxylase activities of the P-450s. Further modeling of the DHEA binding orientation in the substrate-heme pocket of bacterial P-450cam (Iwasaki, M., Darden, T., Pedersen, L., Davis, D. G., Juvonen, R. O., Sueyoshi, T., and Negishi, M. (1993) J. Biol. Chem. 268, 759-762) provides support for the hypothesis that the type of residue at position 117 determines the steroid-substrate specificity of the P-450 depending on the substituent at the C3 position of steroid molecule.

Role of glutamine-61 in the hydrolysis of GTP by p21H-ras: an experimental and theoretical study

The active GTP-bound form of p21ras is converted to the biologically inactive GDP-bound form by enzymatic hydrolysis and this function serves to regulate the wild-type ras protein. The side chain of the amino acid at position 61 may play a key role in this hydrolysis of GTP by p21. Experimental studies that define properties of the Q61E mutant of p21H-ras are presented along with supporting molecular dynamics simulations. We find that under saturating concentrations of GTP the Q61E mutant of p21H-ras has a 20-fold greater rate of intrinsic hydrolysis (kcat = 0.57 min-1) than the wild type. The affinity of the Q61E variant for GTP (Kd = 115 microM) is much lower than that of the wild type. GTPase activating protein does not activate the variant. From molecular dynamics simulations, we find that both the wild type and Q61E mutant have the residue 61 side chain in transient contact with a water molecule that is well-positioned for hydrolytic attack on the gamma phosphate. Thr-35 also is found to form a transient hydrogen bond with this critical water. These elements may define the catalytic complex for hydrolysis of the GTP [Pai et al. (1990) EMBO J. 9, 2351]. Similarly, the G12P mutant, which also has an intrinsic hydrolysis rate similar to the wild type, is found to form the same complex in simulation. In contrast, molecular dynamics analysis of the mutants G12R, G12V, and Q61L, which have much lower intrinsic rates than the wild-type p21, do not show this complex.(ABSTRACT TRUNCATED AT 250 WORDS)

Molecular modeling studies suggest that zinc ions inhibit HIV-1 protease by binding at catalytic aspartates

Human immunodeficiency virus type 1 protease is inhibited in vitro by zinc ions at neutral pH. The binding site of these ions is not known; however, experimental data suggest that binding may occur in the active site. To examine the possibility of zinc binding in the active site, molecular dynamics simulations in the presence and absence of zinc have been carried out to 200 psec. The results are compared with the 2.8-A crystallographic structures of a synthetic HIV-1 protease, and a zinc binding site at the catalytic aspartate residues (Asp-25, Asp-25') is proposed. Molecular dynamics simulations show that the zinc ion remains stably bound in this region, coordinating the carboxylate side chains of both aspartate residues. Interaction with zinc does not disrupt the dimeric structure of the protein or significantly alter the structure of the active site. These data are consistent with experimental studies of HIV-1 protease inhibition by zinc and give strong evidence that this is the binding site that leads to inactivation.

Molecular dynamics simulation of HIV-1 protease in a crystalline environment and in solution

Simulations of the unbound form of the human immunodeficiency virus type 1 protease have been carried out to 200 ps in a crystalline environment and in solution. Solution simulations were performed with and without charge-balancing counterions. The results are compared with the 2.8-A crystallographic structure of Wlodawer et al. [(1989) Science 245, 616], and a proposed model for the solution structure which involves local refolding of the flap regions is presented. The simulations suggest the crystal packing environment of the protease dimer stabilizes the flaps in an extended conformation. Solvation of the dimer leads to local refolding of the flaps which contract toward the active site, forming increased overlap and stronger intersubunit hydrogn bonding at the tips. The degree to which the flaps overlap in solution is observed to depend on the charge state of the system.

Engineering mouse P450coh to a novel corticosterone 15 alpha-hydroxylase and modeling steroid-binding orientation in the substrate pocket

The F209L mutation alters specificity of P450coh from coumarin 7-hydroxylation to 15 alpha-hydroxylation of 11-deoxysteroids such as testosterone and 11-deoxycorticosterone. Neither the wild-type nor F209L exhibits activity toward 11 beta-hydroxysteroids including corticosterone. Mutation of Phe-209 to Asn, however, confers on mutant F209N a high corticosterone 15 alpha-hydroxylase activity. F209V also exhibits low corticosterone 15 alpha-hydroxylase activity; Km and Vmax are 10-fold higher and lower, respectively, than for F209N. The results are consistent with the hypothesis that direct interaction of Asn-209 with 11OH is responsible for high corticosterone 15 alpha-hydroxylase activity. To support this hypothesis, a possible steroid-binding orientation is modeled in the substrate pocket of P450cam. Our weighted homology and constrained alignments map residue 209 of P450coh to Met-184 and Met-191 of P450cam. Energy minimization of corticosterone in the substrate pocket results in the 11OH of the steroid directed toward Met-184 (7 A) and Met-191 (16 A), and in C15 located near the sixth axial position of the heme. The steroid-binding model suggests that the P450cam's substrate pocket may be conserved in the mammalian P450 and can accommodate a steroid molecule, and that residue 209 appears to be located at the critical site that determines the steroid-substrate specificity of a P450 depending on the type of group at the 11-position of steroid molecule.

Simulation of the solution structure of the H-ras p21-GTP complex

An unconstrained simulation of the GTP-bound form of the H-ras protein p21 is performed in an aqueous environment with charge-neutralizing counterions. The simulation is compared to the 1.35-A structure of Pai et al. [(1990) EMBO J. 9, 2351] and a proposed alternate structure, in which the loop at residues 60-65 is modeled into a form which may activate a water molecule for the GTP hydrolysis. The simulation suggests that some protein intermolecular H-bond contacts which are present in the crystal structure are lost in the solvation process and this loss may lead to localized refolding of the molecule. For instance, we find that the gamma-phosphate of the GTP has somewhat weaker contact with the protein in the simulation structure. The antiparallel beta-sheet (residues 38-57) partially melts. The 60-65 loop, which is hypervariable in the X-ray study, is initially relatively distant from the gamma-phosphate region. However, this loop moves so as to sample the space around the gamma-phosphate. For a significant fraction of the simulation time, forms similar to the alternate structure are observed, and a water molecule is localized near the hydrolytic site. The molecular dynamics simulations of p21-GTP in solution support a postulated hydrolysis mechanism for the biological inactivation of the nucleotide complex based on crystallographic data.

Theoretical and experimental measures of DNA helix stability and their relation to sequence specific repair of O6-ethylguanine lesions

Recent work (Breslauer et al. (1986) Proc. Natl. Acad. Sci. (U.S.A.), 83, 3746) has provided a method for calculating empirical thermodynamic quantities for helix to coil transitions from the base sequence of any oligomer. It is shown in this work that the DNA helix binding energy, calculated with the AMBER force field, for 9-mers of the type 5'-GGGXGeYGGG-3', where X and Y are any base and the central Ge is O6-ethylguanine, correlates well with the empirical delta G for helix to strand transitions. The mutation spectrum of ethane methylsulfonate (EMS) in the lacI gene of Escherichia coli can be modeled using the calculated local binding energy but the empirical free energies, enthalpies and melting temperatures predict these levels of repair less well. The relation of the binding energy to the mutation spectrum can be somewhat improved by including entropic effects in a theoretical free energy of binding as given by delta G theoretical identical to delta E binding - T delta S.

A theoretical study of the effect of methylation or ethylation at O6-guanine in the structure and energy of DNA double strands

Quantum and molecular mechanical calculations were employed to examine the effect on binding energies and structure of methylation and ethylation at O6-guanine in double-stranded DNA. Ab initio quantum chemical calculations (STO-3G, 3-21G) were initially used to pseudo-optimize the structure of the 9-methyl derivative of O6-methylguanine. The distal orientation for the O6-methyl group was found to be lower in energy than the proximal orientation. The geometry determined for the distal O6-methyl group was in agreement with recent X-ray work. These results were used in supplementary parameterization of the AMBER molecular mechanics force field necessary for the minimization of DNA double strands containing O6-methylguanosine. Resulting calculations with AMBER on two 5-mer DNA sequences containing the promutagenic G(GM)A subsequence showed that the proximal orientation, while higher in energy in the isolated molecule, is both less disruptive to the DNA double helix and more stable than the distal orientation. Binding energies and degree of destabilization upon methylation were found to be functions of the adjacent bases around a GGA subsequence. Sequence-dependent destabilization could play a role in the repair of alkylated bases. Quantum and molecular mechanics calculations indicate that the O6-methyl and O6-ethylguanines behave energetically in a very similar manner. These calculations suggest that the necessity for the different repair mechanisms for methylation and ethylation lesions cannot be simply explained by energy differences or observed structural differences.

A theoretical study of the binding of polychlorinated biphenyls (PCBs), dibenzodioxins, and dibenzofuran to human plasma prealbumin

Binding energies to human plasma prealbumin using the energy minimization program AMBER are found for a series of polychlorinated biphenyls, dibenzodioxins, and dibenzofuran. Corrections for solvation free energies of the chlorinated analogues lead to estimates of the differential free energies of complex formation. These are compared in a number of cases to known experimental log (KPCB/Kref) values. The theory correctly separates strong, intermediate, and nonbinders. On the basis of calculations, 2,3,7,8-tetrachlorodibenzodioxin and 2,3,7,8-tetrachlorodibenzofuran are predicted to be strong binders, 3,3',5,5'-tetrachlorodiphenoquinone is predicted to be a weak binder, and octachlorodibenzodioxin is predicted to not bind at all. This theoretical model for prealbumin interactions may be of use in estimating the toxic potential of PCBs and related halogenated aromatic hydrocarbons of environmental importance.