Optimized and parallelized implementation of the electronegativity equalization method and the atom-bond electronegativity equalization method.

The most common way to calculate charge distribution in a molecule is ab initio quantum mechanics (QM). Some faster alternatives to QM have also been developed, the so-called "equalization methods" EEM and ABEEM, which are based on DFT. We have implemented and optimized the EEM and ABEEM methods and created the EEM SOLVER and ABEEM SOLVER programs. It has been found that the most time-consuming part of equalization methods is the reduction of the matrix belonging to the equation system generated by the method. Therefore, for both methods this part was replaced by the parallel algorithm WIRS and implemented within the PVM environment. The parallelized versions of the programs EEM SOLVER and ABEEM SOLVER showed promising results, especially on a single computer with several processors (compact PVM). The implemented programs are available through the Web page http://ncbr.chemi.muni.cz/~n19n/eem_abeem.

Molecular dynamics simulations of solvated UDP-glucose in interaction with Mg2+ cations.

Glycosyltransferases are key enzymes involved in biosynthesis of oligosaccharides. Nucleotide-sugars, the glycosyltransferase substrates, serve as activated donors of sugar residues during the enzymatic reaction Although very little is known about the catalytic mechanism of these enzymes, it appears that the catalytic activity in most glycosyltransferases is dependent upon the presence of a divalent cation, for example Mn2+ or Mg2+. It is not known whether the ion is bound to the enzyme before its interaction with the substrate, or if it binds the substrate before the enzymatic reaction to modify its conformation to fit better the active site of the enzyme. We have inspected the latter possibility by running four 2-ns molecular dynamics trajectories on fully solvated UDP-glucose in the presence of Mg2+ ions. Our results indicate that the divalent cation interacts strongly with the nucleotide-sugar in solution, and that it can alter its conformational behavior. It is also shown that a conformation of the pyrophosphate moiety that results in an eclipsed or almost eclipsed orientation of two of the oxygen atoms, and which is found in protein interacting with a nucleotide di- or tri-phosphate X-ray data, is energetically favored. The results are also discussed in light of existing NMR data, and are found to be in a good agreement with them.

Mobility of the active site bound paraoxon and sarin in zinc-phosphotriesterase by molecular dynamics simulation and quantum chemical calculation.

The kinetic data published on phosphotriesterase (PTE), with various complexed metals, clearly indicates that the P=O and P=S bonds of phosphotriester and thiophosphotriester substrates, respectively, are strongly polarized by one or both of the active site complexed metal ions. However, this observation is not consistent with the three-dimensional X-ray crystal structure of zinc-substituted PTE with active site bound substrate analogue diethyl 4-methylbenzylphosphonate. In this structure, the distance between the phosphoryl oxygen and the nearest zinc is 3.4 A, a distance too large to afford strong polarization. In the present paper, the geometry and mobility of various PTE active site-substrate complexes are examined by performing both molecular dynamics (MD) simulations and quantum mechanical calculations. Two known substrates are considered, paraoxon and sarin, although their turnover rates vary about 100-fold. The results indicate that PTE forms a complex with either substrate in which the phosphoryl oxygen becomes strongly coordinated with the less buried zinc atom. It is shown that the geometry of the active site is changed when the protein is immersed in a water bath and relaxed by MD. The most substantial conformational change is the opening of the gateway in a pocket where the location of the leaving group is expected. The opening is observed for the pure enzyme as well as for the enzyme/substrate complexes and it ranges from 11 to 18 A. It is also shown that the pockets, in which the substrate substituents are localized, exhibit different flexibility and interact with the substrate with coordinated conformational adjustments.

Conformational behavior of nucleotide-sugar in solution: molecular dynamics and NMR study of solvated uridine diphosphate-glucose in the presence of monovalent cations.

The nucleotide-sugars are metabolites of primary importance in the biosynthesis of polysaccharides and glycoconjugates since they serve as sugar donors in the reactions of glycosyltransferases, enzymes that displays a high specificity for both donors and acceptors. In order to determine the conformational behavior of uridinediphosphoglucose in dilute aqueous solution that includes a physiologically relevant concentration of salt, parallel NMR and molecular modeling investigations have been conducted. Nine molecular dynamics trajectories of 3 ns each were calculated in presence of explicit water and monovalent cations with the use of the AMBER force field with recently developed energy parameters for nucleotide-sugars (P. Petrova, J. Koca, and A. Imberty, Journal of American Chemical Society, 1999, vol. 121, pp. 5535-5547). Theoretical nuclear Overhauser effect data were calculated from these simulations using a model-free approach that takes into account internal motions. Comparison of theoretical and experimental data gives excellent agreement for the region surrounding the glucose-phosphate linkage including the pyrophosphate linkage itself. Less satisfactory agreement is obtained for the ribose ring and the base orientations. On the whole, both NMR and molecular dynamics simulations predict the molecule to be flexible, and to visit a large number of conformations while maintaining an extended overall shape.

TRITON: graphic software for rational engineering of enzymes.

Engineering of the catalytic properties of enzymes requires knowledge about amino acid residues interacting with the transition state of the substrate. TRITON is a graphic software package for modelling enzymatic reactions for the analysis of essential interactions between the enzyme and its substrate and for in silico construction of protein mutants. The reactions are modelled using semi-empirical quantum-mechanic methods and the protein mutants are constructed by homology modelling. The users are guided through the calculation and data analysis by wizards.

A-like guanine-guanine stacking in the aqueous DNA duplex of d(GGGGCCCC).

We have used CD spectroscopy, NMR spectroscopy and unrestrained molecular dynamics to study conformational properties of a DNA duplex formed by the self-complementary octamer d(GGGGCCCC). Its unusual CD spectrum contains features indicating A-like stacking of half of the bases, whereas the other half stack in a B-like fashion. Unrestrained molecular dynamics simulations converged to a stable B-like double-helix of d(GGGGCCCC). However, the double-helix contained a central hole whose size was half of that occurring in structure A. In the canonical structure B, the hole does not exist at all because the base-pairs cross the double-helix centre. The cytosine bases were stacked in the duplex of d(GGGGCCCC) as in structure B, while stacking of the guanine bases displayed features characteristic for structure A. NMR spectroscopy revealed that the A-like guanine-guanine stacking was accompanied by an increased tendency of the deoxyribose rings attached to the guanine bases to be puckered in an A-like fashion. Otherwise, the duplex of d(GGGGCCCC) showed no clash, no bend and no other significant deviation from structure B. The present analysis demonstrates a remarkable propensity of the guanine runs to stack in an A-like fashion even within the B-DNA framework. This property explains why the oligo(dG). oligo(dC) tracts switch into structure A so easily. Secondly, this property may influence replication, because structure A is replicated more faithfully than structure B. Thirdly, the oligo(dG) runs might have played an important role in early evolution, when DNA took on functions that originally evolved on RNA. Fourthly, the present study extends the vocabulary of DNA secondary structures by the heteronomous duplex of d(GGGGCCCC) in which the B-like strand of oligo(dC) is bound to the A-like strand of oligo(dG).

Molecular dynamics of DNA quadruplex molecules containing inosine, 6-thioguanine and 6-thiopurine.

The ability of the four-stranded guanine (G)-DNA motif to incorporate nonstandard guanine analogue bases 6-oxopurine (inosine, I), 6-thioguanine (tG), and 6-thiopurine (tI) has been investigated using large-scale molecular dynamics simulations. The simulations suggest that a G-DNA stem can incorporate inosines without any marked effect on its structure and dynamics. The all-inosine quadruplex stem d(IIII)(4) shows identical dynamical properties as d(GGGG)(4) on the nanosecond time scale, with both molecular assemblies being stabilized by monovalent cations residing in the channel of the stem. However, simulations carried out in the absence of these cations show dramatic differences in the behavior of d(GGGG)(4) and d(IIII)(4). Whereas vacant d(GGGG)(4) shows large fluctuations but does not disintegrate, vacant d(IIII)(4) is completely disrupted within the first nanosecond. This is a consequence of the lack of the H-bonds involving the N2 amino group that is not present in inosine. This indicates that formation of the inosine quadruplex could involve entirely different intermediate structures than formation of the guanosine quadruplex, and early association of cations in this process appears to be inevitable. In the simulations, the incorporation of 6-thioguanine and 6-thiopurine sharply destabilizes four-stranded G-DNA structures, in close agreement with experimental data. The main reason is the size of the thiogroup leading to considerable steric conflicts and expelling the cations out of the channel of the quadruplex stem. The G-DNA stem can accommodate a single thioguanine base with minor perturbations. Incorporation of a thioguanine quartet layer is associated with a large destabilization of the G-DNA stem whereas the all-thioguanine quadruplex immediately collapses.

Stability of complexes of aromatic amides with bromide anion: quantitative structure--property relationships

Most of the theoretical studies published to-date on the structural and electronic properties of supramolecules have been devoted to the neutral or cationic complexes, while little is known about anionic systems. A detailed theoretical study of the interaction between simple aromatic amides and the bromide anion has recently been published (Cajan, M.; Stibor, I.; Koca, J. J. Phys. Chem. A 1999, 103, 3778). The present work focuses on the structural and physicochemical parameters of simple aromatic amides related to their ability to form the 1:1 complex with a bromide anion. A quantitative structure-property relationships (QSPR) model for the prediction of association constants is proposed. The model based on 22 complexes and nine molecular descriptors explained 96% (84% cross-validated) of the variance in association constants. The descriptors employed in this model included parameters for the characterization of conformational behavior and the 3D structure of amide molecules, distribution of electron density on the amidic functional group, and parameters for substitution on aromatic units. The quantitative structure-property relationship approach predicted the association constants with comparable quality, but significantly lower computational demand, than molecular modeling or standard quantum chemistry calculations.

Single-coordinate-driving method for molecular docking: application to modeling of guest inclusion in cyclodextrin.

An extension of the computer program CICADA has been developed that allows us to use the single-coordinate-driving (SCD) method for flexible molecular docking. The docking procedure is composed of three independent space rotations, three independent translations, and the torsions selected by the user. One of the coordinates is driven; the other coordinates are relaxed. This procedure follows low-energy wells on the potential energy surface of the entire system. The program allows us to dock more than one ligand molecule to the receptor. We ran two test examples, docking N,N-dimethylformamide into alpha-cyclodextrin and R-phenoxypropionic acid into beta-cyclodextrin. The test examples showed that the SCD approach is able to overcome high-energy barriers and to cover the entire box within which the search is performed. The limitations of molecular dynamics docking in comparison with our approach also are discussed. The philosophy of the newly developed approach is not only to find the best dock for the receptor-ligand(s) system, but also to describe all the important binding modes and provide a good starting point for studying the dynamics within the cavity during the docking process.

Docking-based development of purine-like inhibitors of cyclin-dependent kinase-2.

The cell division cycle is controlled by cyclin-dependent kinases (cdk), which consist of a catalytic subunit (cdk1-cdk8) and a regulatory subunit (cyclin A-H). Purine-like inhibitors of cyclin-dependent kinases have recently been found to be of potential use as anticancer drugs. Rigid and flexible docking techniques were used for analysis of binding mode and design of new inhibitors. X-ray structures of three (ATP, olomoucine, roscovitine) cdk2 complexes were available at the beginning of the study and were used to optimize the docking parameters. The new potential inhibitors were then docked into the cdk2 enzyme, and the enzyme/inhibitor interaction energies were calculated and tested against the assayed activities of cdk1 (37 compounds) and cdk2 (9 compounds). A significant rank correlation between the activity and the rigid docking interaction energy has been found. This implies that (i) the rigid docking can be used as a tool for qualitative prediction of activity and (ii) values obtained by the rigid docking technique into the cdk2 active site can also be used for the prediction of cdk1 activity. While the resulting geometries obtained by the rigid docking are in good agreement with the X-ray data, the flexible docking did not always produce the same inhibitor conformation as that found in the crystal.

A conformational study of the xyloglucan oligomer, XXXG, by NMR spectroscopy and molecular modeling.

A structural study of the XXXG xyloglucan heptasaccharide (X = alpha-D-Xylp(1 --> 6)-beta-D-Glcp and G = beta-D-Glcp) isolated from apple fruit has been undertaken with nmr and molecular mechanics methods. Quantitative 400 MHz nmr data including nuclear Overhauser effect spectroscopy (NOESY) volumes were recorded at both 6 and 20 degrees C. In spite of severe overlapping of resonances, it was possible to estimate summed NOEs for the majority of the anomeric and glucosyl methylene protons. An ensemble-average population of preferred geometries has been established with the CICADA conformational searching algorithm associated with the MM3 force field. Comparison of the theoretical data obtained by back-calculation of the NOESY volumes from the ensemble-average distance matrix program and motional models based on the Stokes-Einstein-Debye relation satisfactorily reproduce the experimental data. Conformational averaging about the mainchain glycosidic linkages includes both the syn and anti conformers and a minor gauche-gauche population is highly probable. The theoretical data overestimate the syn preference of the Glc(c) --> Glc(b) linkage as well as the Glc(c) GT rotamer population. Finally, both the motional models and the conformational search indicate a fairly rigid backbone and greater flexiblity for the xylose side chains.

An A-type double helix of DNA having B-type puckering of the deoxyribose rings.

DNA usually adopts structure B in aqueous solution, while structure A is preferred in mixtures of trifluoroethanol (TFE) with water. However, the octamer d(CCCCGGGG) and other d(C(n)G(n)) fragments of DNA provide CD spectra that suggest that the base-pairs are stacked in an A-like fashion even in aqueous solution. Yet, d(CCCCGGGG) undergoes a cooperative TFE-induced transition into structure A, indicating that an important part of the aqueous duplex retains structure B. NMR spectroscopy shows that puckering of the deoxyribose rings is of the B-type. Hence, combination of the information provided by CD spectroscopy and NMR spectroscopy suggests an unprecedented double helix of DNA in which A-like base stacking is combined with B-type puckering of the deoxyribose rings. In order to determine whether this combination is possible, we used molecular dynamics to simulate the duplex of d(CCCCGGGG). Remarkably, the simulations, completely unrestrained by the experimental data, provided a very stable double helix of DNA, exhibiting just the intermediate B/A features described above. The double helix contained well-stacked guanine bases but almost unstacked cytosine bases. This generated a hole in the double helix center, which is a property characteristic for A-DNA, but absent from B-DNA. The minor groove was narrow at the double helix ends but wide at the central CG step where the Watson-Crick base-pairs were buckled in opposite directions. The base-pairs stacked tightly at the ends but stacking was loose in the duplex center. The present double helix, in which A-like base stacking is combined with B-type sugar puckering, is relevant to replication and transcription because both of these phenomena involve a local B-to-A transition.

TRITON: in silico construction of protein mutants and prediction of their activities.

MOTIVATION: One of the objectives of protein engineering is to propose and construct modified proteins with improved activity for the substrate of interest. Systematic computational investigation of many protein variants requires the preparation and handling of a large number of data files. The type of the data generated during the modelling of protein variants and the estimation of their activities offers the possibility of process automatization. RESULTS: The graphical program TRITON has been developed for modelling protein mutants and assessment of their activities. Protein mutants are modelled from the wild type structure by homology modelling using the external program MODELLER. Chemical reactions taking place in the mutants active site are modelled using the semi-empirical quantum mechanic program MOPAC. Semi-quantitative predictions of mutants activities can be achieved by evaluating the changes in energies of the system and partial atomic charges of active site residues during the reaction. The program TRITON offers graphical tools for the preparation of the input data files, for calculation and for the analysis of the generated output data. AVAILABILITY: The program TRITON can run under operating systems IRIX, Linux and NetBSD. The software is available at http://www.chemi.muni.cz/lbsd/triton.ht ml.

Analysis of the reaction mechanism and substrate specificity of haloalkane dehalogenases by sequential and structural comparisons.

Haloalkane dehalogenases catalyse environmentally important dehalogenation reactions. These microbial enzymes represent objects of interest for protein engineering studies, attempting to improve their catalytic efficiency or broaden their substrate specificity towards environmental pollutants. This paper presents the results of a comparative study of haloalkane dehalogenases originating from different organisms. Protein sequences and the models of tertiary structures of haloalkane dehalogenases were compared to investigate the protein fold, reaction mechanism and substrate specificity of these enzymes. Haloalkane dehalogenases contain the structural motifs of alpha/beta-hydrolases and epoxidases within their sequences. They contain a catalytic triad with two different topological arrangements. The presence of a structurally conserved oxyanion hole suggests the two-step reaction mechanism previously described for haloalkane dehalogenase from Xanthobacter autotrophicus GJ10. The differences in substrate specificity of haloalkane dehalogenases originating from different species might be related to the size and geometry of an active site and its entrance and the efficiency of the transition state and halide ion stabilization by active site residues. Structurally conserved motifs identified within the sequences can be used for the design of specific primers for the experimental screening of haloalkane dehalogenases. Those amino acids which were predicted to be functionally important represent possible targets for future site-directed mutagenesis experiments.

Differences in conformational behavior of ATA and TAT sequences in single strand DNA trimer.

The conformational behavior of single strand (ss) TAT and ATA trimers of DNA have been studied by computational chemistry tools including CICADA software interfaced with AMBER molecular mechanics and dynamics. The Single-Coordinate-Driving (SCD) method has been used in conjunction with molecular dynamics simulated annealing. It has been revealed that the conformational flexibility of each sequence differs substantially from the other one. Four common conformational families have been found for both trimers. These are: helical, reverse-stacked (base 3), half-stacked (base 3), reverse-stacked (base 1). However, the energies of conformers representing the families are different for both the studied systems. An additional conformational family, bulged, has been found for ss(ATA), while ss(TAT) has been found also in half-stacked (base 1) conformation. In general, ss(TAT) exhibits a higher number of low energy conformations while ss(ATA) shows one interesting low energy conformational interconversion between reverse-stacked (A3) family and half-stacked (A3) family. The high conformational variability of the trimers has been confirmed by flexibility analysis and by molecular dynamics simulations, which have also shown the conformational stability of single conformational families. It has been concluded that the methodology used is able to provide a very detailed picture of the conformational space of these molecules.

Computational site-directed mutagenesis of haloalkane dehalogenase in position 172.

The application of molecular modelling and quantum-chemistry calculations for the 'computational site-directed mutagenesis' of haloalkane dehalogenase is described here. The exhaustive set of single point mutants of haloalkane dehalogenase in position 172 was constructed by homology modelling. The ability of substituting residues to stabilize the halide ion formed during the dehalogenation reaction in the enzyme active site was probed by quantum-chemical calculations. A simplified modelling procedure was adopted to obtain informative results on the potential activity of mutant proteins in a sufficiently short period of time, which, in the future, could be applicable for making bona fide predictions of mutants' activity prior to their preparation in the laboratory. The reaction pathways for the carbon-halide bond cleavage were calculated using microscopic models of wild type and mutant proteins. The theoretical parameters derived from the calculation, i.e. relative energies and selected atomic charges of educt, product and transition state structures, were statistically correlated with experimentally determined activities. The charge difference of educt and product on the halide-stabilizing hydrogen atom of residue 172 was the best parameter to distinguish protein variants with high activity from mutant proteins displaying a low activity. All mutants with significant activity in the experiment were found to have this parameter one order of magnitude higher than mutants with low activity. The results obtained are discussed in the light of the practical application of this methodology for the prediction of potentially active protein variants. Further automation of the modelling procedure is suggested for combinatorial screening of the large number of protein variants. Coupling of the dehalogenation reaction with hydrogenation of the halide ion formed during the reaction in the enzyme active site was proposed as a possible way to improve the catalytic activity of the haloalkane dehalogenase of Xanthobacter autotrophicus GJ10.

Travelling through conformational space: an approach for analyzing the conformational behaviour of flexible molecules.

The applications of the single-co-ordinate-driving (SCD) method in conformational analysis of flexible molecules have been discussed. SCD can best be characterised as travelling through low energy areas of the conformational space. It has been shown that SCD provides detailed information about the conformational behaviour of small and middle sized flexible molecules. It has been demonstrated that SCD may fail for molecules which are more rigid but still conformationally interesting, for example, AUG trimer of RNA. It has been found out that the search problems are eliminated when SCD is coupled with simulated annealing (SCD-SA). Both SCD and SCD-SA methods may be recognised as successful tools for analysis of conformational space. It has been demonstrated, how SCD and SCD-SA results can be used as a background to discover and analyse correlated conformational processes, to quantify molecular conformational flexibility, and to provide an appropriate background for an efficient free energy simulation.

Repositioning the catalytic triad aspartic acid of haloalkane dehalogenase: effects on stability, kinetics, and structure.

Haloalkane dehalogenase (DhlA) catalyzes the hydrolysis of haloalkanes via an alkyl-enzyme intermediate. The covalent intermediate, which is formed by nucleophilic substitution with Asp124, is hydrolyzed by a water molecule that is activated by His289. The role of Asp260, which is the third member of the catalytic triad, was studied by site-directed mutagenesis. Mutation of Asp260 to asparagine resulted in a catalytically inactive D260N mutant, which demonstrates that the triad acid Asp260 is essential for dehalogenase activity. Furthermore, Asp260 has an important structural role, since the D260N enzyme accumulated mainly in inclusion bodies during expression, and neither substrate nor product could bind in the active-site cavity. Activity for brominated substrates was restored to D260N by replacing Asn148 with an aspartic or glutamic acid. Both double mutants D260N+N148D and D260N+N148E had a 10-fold reduced kcat and 40-fold higher Km values for 1,2-dibromoethane compared to the wild-type enzyme. Pre-steady-state kinetic analysis of the D260N+N148E double mutant showed that the decrease in kcat was mainly caused by a 220-fold reduction of the rate of carbon-bromine bond cleavage and a 10-fold decrease in the rate of hydrolysis of the alkyl-enzyme intermediate. On the other hand, bromide was released 12-fold faster and via a different pathway than in the wild-type enzyme. Molecular modeling of the mutant showed that Glu148 indeed could take over the interaction with His289 and that there was a change in charge distribution in the tunnel region that connects the active site with the solvent. On the basis of primary structure similarity between DhlA and other alpha/beta-hydrolase fold dehalogenases, we propose that a conserved acidic residue at the equivalent position of Asn148 in DhlA is the third catalytic triad residue in the latter enzymes.

CICADA interface with AMBER. An application on oligonucleotides and their fragments.

The potential energy hypersurfaces (PES) of several nucleotide fragments were analyzed by the conformational search algorithm CICADA interfaced with the molecular mechanics program AMBER, version 4.0. The low energy conformers for dimethylphosphate, dinucleoside monophosphate fragments, and deoxyadenosine are described. Calculated relative and absolute flexibilities of single conformers, molecular fragments as well as entire molecules are introduced. The comparison of the results with the literature data shows good ability of the CICADA-AMBER combination to describe conformational space. It is revealed that the number of low energy conformers as well as flexibility decreases as the size of the molecule increases. The conformational behavior of freely rotatable single bonds, specially those within a sugar ring, is more "sharp" in larger structures.

Travelling on the potential energy surfaces of carbohydrates: comparative application of an exhaustive systematic conformational search with an heuristic search.

The calculated ensembles found by a heuristic conformational search algorithm, CICADA, for three small carbohydrates, ethyl beta-lactoside, methyl alpha-D-galactoside, and methyl beta-D-galactoside, are evaluated in terms of their ability to reproduce time-averaged optical rotation and NMR data. A unique dynamic model for methyl beta-D-galactoside has been obtained by fitting experimental NOESY volumes to the theoretical ones elaborated from the CICADA ensemble internuclear distances with the model-free formalism. In the case of ethyl beta-lactoside, the CICADA ensemble is compared to that of an exhaustive systematic grid-search method. The CICADA algorithm proved to be a very efficient method to find most of the important minima on even very complex potential energy surfaces, and the spectral quality of the CICADA ensemble was found to be of equal quality, if not superior, to that of the exhaustive systematic grid-search method. The CICADA algorithm has several advantages over other conformational search algorithms: (1) It has polynomial dependence of dimensions on computer time in contrast to the grid search, which has exponential dependence, (2) the conformations found are free of artificial harmonic constraint potentials, (3) it passes all barriers amongst families of conformations on conformational hypersurface but spends almost all its time in the essential highly populated areas, (4) the inherent properties of the algorithm make rigorous minimization criteria superfluous and provide good convergence behavior, and (5) as an important spin-off, it provides low-energy interconversion pathways that can, amongst others, be used for estimating adiabatic rotational barriers.