1471-2105-9-438, bioinformatyka, artykuly

//-->.pos {position:absolute; z-index: 0; left: 0px; top: 0px;}BMC BioinformaticsMethodology articleBioMedCentralOpen AccessAMMOS: Automated Molecular Mechanics Optimization tool forinsilicoScreeningTania Pencheva1,2, David Lagorce1, Ilza Pajeva2, Bruno O Villoutreix1andMaria A Miteva*1Address:1INSERM U648, Bioinformatics-MTI University Paris Diderot, 5 rue Marie-Andrée Lagroua, 75205 Paris Cedex 13, France and2Centre ofBiomedical Engineering "Prof. Ivan Daskalov", Bulgarian Academy of Sciences, 105, Akad. Georgi Bonchev Str., 1113 Sofia, BulgariaEmail: Tania Pencheva - tania.pencheva@clbme.bas.bg; David Lagorce - david.lagorce@univ-paris-diderot.fr; Ilza Pajeva - pajeva@bio.bas.bg;Bruno O Villoutreix - bruno.villoutreix@univ-paris-diderot.fr; Maria A Miteva* - maria.miteva@univ-paris-diderot.fr* Corresponding authorPublished: 16 October 2008BMC Bioinformatics2008,9:438doi:10.1186/1471-2105-9-438Received: 18 April 2008Accepted: 16 October 2008This article is available from: http://www.biomedcentral.com/1471-2105/9/438© 2008 Pencheva et al; licensee BioMed Central Ltd.This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0),which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.AbstractBackground:Virtual orin silicoligand screening combined with other computational methods isone of the most promising methods to search for new lead compounds, thereby greatly assistingthe drug discovery process. Despite considerable progresses made in virtual screeningmethodologies, available computer programs do not easily address problems such as: structuraloptimization of compounds in a screening library, receptor flexibility/induced-fit, and accurateprediction of protein-ligand interactions. It has been shown that structural optimization of chemicalcompounds and that post-docking optimization in multi-step structure-based virtual screeningapproaches help to further improve the overall efficiency of the methods. To address some of thesepoints, we developed the program AMMOS for refining both, the 3D structures of the smallmolecules present in chemical libraries and the predicted receptor-ligand complexes throughallowing partial to full atom flexibility through molecular mechanics optimization.Results:The program AMMOS carries out an automatic procedure that allows for the structuralrefinement of compound collections and energy minimization of protein-ligand complexes using theopen source program AMMP. The performance of our package was evaluated by comparing thestructures of small chemical entities minimized by AMMOS with those minimized with the Triposand MMFF94s force fields. Next, AMMOS was used for full flexible minimization of protein-ligandscomplexes obtained from a mutli-step virtual screening. Enrichment studies of the selected pre-docked complexes containing 60% of the initially added inhibitors were carried out with or withoutfinal AMMOS minimization on two protein targets having different binding pocket properties.AMMOS was able to improve the enrichment after the pre-docking stage with 40 to 60% of theinitially added active compounds found in the top 3% to 5% of the entire compound collection.Conclusion:The open source AMMOS program can be helpful in a broad range ofin silicodrugdesign studies such as optimization of small molecules or energy minimization of pre-dockedprotein-ligand complexes. Our enrichment study suggests that AMMOS, designed to minimize alarge number of ligands pre-docked in a protein target, can successfully be applied in a final post-processing step and that it can take into account some receptor flexibility within the binding sitearea.Page 1 of 15(page number not for citation purposes)BMC Bioinformatics2008,9:438http://www.biomedcentral.com/1471-2105/9/438BackgroundStructure-based virtual ligand screening (SBVLS) allows toinvestigate thousands or millions of molecules against abiomolecular target [1,2], and as such it plays an increas-ingly important role in modern drug discovery programs.For example, numerous SBVLS methods employing dock-ing and scoring have been developed to assist the discov-ery of hit compounds and their optimization to leads [3-5]. These methods orient and score small molecules in aprotein-binding site, searching for shape and chemicalcomplementarities. Many novel active compounds actingon key therapeutic targets have been found through com-bining SBVLS and in vitro screening experiments [5,6].Despite the considerable progresses achieved these recentyears, several problems are still present in most of the cur-rently available SBVLS packages. Among the most criticalis the flexibility of the receptors that frequently changetheir conformations upon ligand binding. Several meth-ods have been developed to attempt to take into consider-ation receptor flexibility during docking/scoring [2,7-10],however, this is still very challenging because the numberof conformations rises exponentially with the number ofrotatable bonds and the full sampling of all possible con-formations is not feasible for a large number of protein-ligand complexes.Further the correct prediction of receptor-ligand bindingenergies [11,12] and accurate ranking of the compoundswith respect to their estimated affinities to a targetremains highly challenging. Thus it is still difficult to dis-criminate bioactive compounds from false positives[13,14] despite recent efforts to improve enrichment via,for instance, docking on different protein targets [15] orthrough optimized or new scoring functions [12,16,17].In addition, and among the many players that are impor-tant in SBVLS computations, the quality of the screenedchemical libraries has also been shown to be important inorder to correctly predict the bound ligand-conformationsand for ranking [18,19]. Within this context, furtherrefinements and optimization of VLS docking-scoringmethods are needed.Recently it has been suggested that post-docking optimi-zation, either after conventional docking-scoring proce-dures or after hierarchical VLS protocols [20-23] may helpto further improve both, the docking pose and the scor-ing, and as such the overall efficiency of SBVLS experi-ments. Recent examples of docked poses and enrichmentimprovements after post-docking energy minimizationsupport this view [19,24-27].In the present study, we propose a new open source pro-gram, named AMMOS, which addresses some of the pre-and post-processing problems associated with SBVLScomputations, through molecular mechanics (MM) mod-eling. AMMOS executes an automatic procedure for: (1)energy minimization of pre-docked protein-ligand com-plexes allowing partial or full atom flexibility from both,the ligand and the receptor sides and (2) structural opti-mization of chemical compounds present in the screeninglibraries prior to docking experiments. MM is currently avery reliable approach to model protein-receptor interac-tions in a physically realistic manner [26-28] since it canaccount for local flexibility adjustments from both, theprotein and the ligand although conformational explora-tion is not possible if large conformational changes occur.It is indeed reasonable to apply such framework instead ofmore computer demanding simulations (for instancemolecular dynamics) in large-scale applications involvingthe handling of thousands of compounds. In conven-tional MM studies, the bonded interactions include thebonds, bond angles and dihedral terms while the non-bonded interactions involve the van der Waals term repre-sented by the Lennard-Jones (LJ) 6–12 potential, and elec-trostatic interactions, often treated by Coulombicpotential computed between point charges centered onrelevant atoms. AMMOS proposes relatively fast energyminimization by making use of the full-featured molecu-lar mechanics program AMMP [28-31]. AMMP is availableupon GNU license and has been recently implemented inthe well-known OpenGL molecular modeling packageVEGA [32]. In particular, VEGA implements AMMP forenergy minimization with all the available optimizationmethods. However, to the best of our knowledge, VEGAcan not minimize chemical libraries nor a large number ofpre-docked protein-ligand complexes. AMMP has severaladvantages that make it relevant for the present applica-tions such as: a fast multipole method for including allatoms in the calculation of long-range potentials androbust structural optimizers. AMMP has a flexible choiceof several bonded and non-bonded potentials and per-mits analysis of individual energy terms. An additionaladvantage of AMMP is that it allows straightforward intro-duction of non-standard polymer linkages and non-standard amino-acid residues as well as manipulation ofboth small molecules and macromolecules including pro-teins, nucleic acids and other polymers. Furthermore,extensive benchmarking of AMMP has been performedhighlighting its accuracy in term of energy minimizationof proteins (i.e., the changes of atomic positions afterminimization have been shown to be within the range ofexperimentally obtained variations in different crystalforms [33] and the calculated protein-ligand interactionenergies have been successfully correlated with measuredligand binding affinities [34]). Up to now, AMMP hasbeen successfully applied in numerous modeling studiesof proteins and protein-ligand complexes [28,29,34,35]but it has not been used thus far for SBVLS computationsas numerous implementations would be necessary toapply it automatically on thousands of small molecules,Page 2 of 15(page number not for citation purposes)BMC Bioinformatics2008,9:438http://www.biomedcentral.com/1471-2105/9/438either present in a database or docked in a receptor bind-ing site.In this article we describe the development, implementa-tion and evaluation of the AMMOS approach for auto-matic energy minimization of protein-ligand complexesor of small organic molecules. First energy minimizationwith AMMOS was validated by refinement of a largechemical library and by comparison of small moleculeoptimizations with two well established force fields avail-able in the package SYBYL [36]. The efficiency of AMMOSfor optimization of pre-docked protein-ligand complexeswas then examined through performing calculations ontwo protein targets, namely estrogen receptor (ER) andneuraminidase (NA). These two proteins were selectedbecause they are often used in SBVLS studies and becausetheir binding pockets display rather different physico-chemical properties and geometries. We used a multi-stepSBVLS protocol to generate protein-ligand complexes.Our first VLS stage employed a search for satisfactoryshape complementarity allowing rational reduction of thesize of the initial chemical library. The pre-docked com-plexes subjected to AMMOS minimization were generatedvia the second step involving flexible-ligand docking.Finally we tested AMMOS as a final re-scoring engine onthe selected pre-docked protein-ligand complexes forboth, ER and NA. The obtained results are promising andsuggest that AMMOS can be successfully applied as a post-processing tool to improve enrichment.correlation with experimental dipole moments. MOPSAparameter set has been merged with the AMMP parametersetsp4to generate the new standard force field setsp5[39].We tested several minimization methods as implementedin AMMP. The steepest descent method [40] (first order)that uses the first derivatives of the energy function to finda local minimum, is available in AMMP, but this does notnecessarily produce the fastest convergence. The Poliak-Ribeire Conjugate Gradient method (first order) [40], per-forming a search along conjugate directions, can producegenerally faster convergence, is also implemented inAMMP. For the second-order minimization methods, theBroyden-Fletcher-Goldfarb-Shanno (BFGS) approach[41] that belongs to the Quasi-Newton methods [40], isavailable in AMMP. The non-derivative polytope simplexmethod [41], as well as a genetic algorithm (GA) [42](evolutionary algorithms that perform directed randomsearch to find the optimal solution in a complex multidi-mensional space) are also implemented in AMMP.Energy minimization protocolsEnergy minimization of small organic molecules using AMMOSWe explored both AMMP force fields,sp4andsp5,on a setof four small molecules structures generated by OMEGA2.0 [43] (see Compound collections below). To speed-upligand parameterization, partial charges on ligand atomswere assigned with the Gasteiger-Marsili method [44]using the OpenBabel package [45]. The maximumnumber of iterations was set to 5000. All calculations wereperformed with a convergence value set to 0.01 or 0.02kcal.mol-1.Å-1, and no essential differences were observed.Thus, we chose 0.02 kcal.mol-1.Å-1as convergence crite-rion to reduce the computational time. The number ofiterations required to reach convergence with the conju-gate gradient method (our results demonstrate that thisapproach is the most efficient, see in the Results section)varied between 300 and 1600 for different small com-pounds. Therefore for further analysis, two protocols wereassessed: one protocol employing two subsequent steps of500 iterations and one with 1000 iterations. The mini-mized and initial structures were compared based on theRMSD values between the non-hydrogen atoms using theSuperimpose option of the InsightII molecular modelingpackage [46]. We should notice that small molecule min-imization results depend on initial conformations andthis holds for any MM minimization engine. Problemswith possible biases due to the starting conformation andexploration of other conformers can be circumvented bythe use of our multiconformer generator Multiconf-DOCK [47] or of OMEGA prior to minimizing small mol-ecules.MethodsAMMP Molecular MechanicsAMMOS makes use of the molecular simulation packageAMMP [30], a program that can easily be embedded inother packages. AMMP allows to introduce standard ornon-standard polymer linkages (ensured by the programPREAMMP included in the package AMMP), unusual lig-ands or non-standard residues, as well as to complete par-tial protein structures. AMMP incorporates a fastmultipole algorithm for the efficient calculation of long-range forces thereby allowing evaluation of non-bondedterms without the use of a cutoff radius and increasing thespeed, making calculations comparable to a standardtreatment with a 8–10 Å radius cutoff [28]. The AMMPforce field [28] is developed on the basis of the UFF poten-tial set [37] and the AMBER partial charges [38]. The ini-tial UFF set has been optimized for biological molecules[28,29] in order to improve the agreement with experi-mental geometry and spectral data. The first developedAMMP force field was the setsp4[28]. Lately, Bagossi andco-authors [39] proposed theModified Parameter SetforAMMP (MOPSA) with improved generation of partialcharges for a wider range of compounds especiallyadapted for modeling of macromolecules, where the elec-trostatic parameters have been modified to achieve betterPage 3 of 15(page number not for citation purposes)BMC Bioinformatics2008,9:438http://www.biomedcentral.com/1471-2105/9/438Energy minimization of small organic molecules using other forcefieldsTwo MM minimization methods implemented in the pro-gram SYBYL were applied on the same set of small mole-cules in order to compare the minimization resultsobtained with AMMOS. The initial structures of the foursmall compounds, generated by OMEGA (see for detailsCompound collections), were optimized by two forcefields: the Tripos force field (Tff) [36] and MMFF94s [48].Tff minimization was performed with Gasteiger-Huckelcharges, and MMFF94s, with MMFF charges. For bothforce fields, the following settings were used: distance-dependent dielectric function; non-bonded cutoff 8.0 Å;0.02 kcal.mol-1.Å-1energy gradient convergence criterion;simplex initial optimization. For comparison, several runswere performed with 0.01 kcal.mol-1.Å-1convergence andno essential differences in the resulting geometries wererecorded. Two gradient methods were experimented:Powell and conjugate gradient. Powell method [49]belongs to the conjugate gradient family of minimizationmethods. It is also more tolerant to inexact line searches.As a result, it is faster than the conjugate gradient methodand is well-suited for a wide variety of problems [49]. Thenumber of iterations was set to 5000 for both methods, inall cases however, the convergence was reached wellbelow this number. Because the Powell and conjugate gra-dient results were quite similar, we report here only thedata obtained with the Powell method. The minimizedand initial structures were compared based on the RMSDvalues between the non-hydrogen atoms using the Matchoption in SYBYL. The calculations with SYBYL were per-formed on a Silicon Graphics Octane 2 (R12000) runningunder IRIX 6.5.Docking and scoring protocolThe pre-docked protein-ligand complexes subjected toenergy minimization with AMMOS were generated via amulti-step docking-scoring protocol with DOCK6 [50].DOCK6 accomplishes a sphere-matching algorithm to fitligand atoms to spheres representing a negative image ofthe receptor-binding site. We used the program DMS [51]to compute the molecular surface of the receptor. Theoverlapping spheres within a radius of 4 Å were generatedon the protein binding site surface using the programSPHGEN [52]. Sphere clusters within 6 Å to a referenceligand were retained for ER and 4 Å for NA (i.e., NA pos-sesses a very open and flat binding site and it is importantto limit the search on the active site area for this kind oftarget). The first docking step was carried out using rigidbody-docking with DOCK6 applying the MS-DOCK pro-tocol [47] over a compound collection of 37970 mole-cules (ADME/Tox filtered ChemBridge Diversity set)present in a multi-conformer state without initial minimi-zation with AMMOS (see Compound collections). Thisstage assesses only shape complementarity and therefore,multi-conformer structures for the small compounds areneeded in order to perform this fast "geometric" filteringstep. For the positioning of the ligand in the binding site,we applied the faster manual match (see [50] and a maxi-mum of 500 orientations). As mentioned above, in ourcalculations, the scores measured only the steric comple-mentarity by use of the contact scoring function thatcounts the number of receptor-ligand contacts within a4.5 Å distance from the ligand atoms. Each clash penal-ized the score by 30. The allowed bump overlaps werechosen to be 0.75 for NA and 0.50 for ER. These valueswere selected according to our previous observations [47].We have seen that for large and flat cavities binding rela-tively small ligands like in the case of NA, a bump overlapof 0.75 improves the enrichment after a rigid docking pro-cedure. On the contrary, when large ligands fill well thebinding site, a bump overlap of 0.50 is preferable, a situa-tion encountered with ER (see for details [47]).Secondly, the retrieved non-minimized top ranked com-pounds (30–50 % of the library containing at least 60% ofthe actives) were directly re-docked using a flexible dock-ing mode (i.e., flexibility from the ligand side) imple-mented in DOCK6 and employing the incremental builtalgorithm "anchor-first" [53] with our optimized parame-ters to better handle ligand flexibility [47]. We used a max-imum of 1000 orientations for the anchor fragment. Tospeed-up the calculations, we set 50 configurations percycle for the growth of the ligands. We applied 20 simplexminimization steps to each growth step. All docked mole-cules were ranked using the standard DOCK score involv-ing soft van der Waals and distance-dependentelectrostatic potentials. Finally for each ligand we savedup to 10 best scored conformers with a RMSD of 0.8 Å forsubsequent minimization with AMMOS.Dataset preparationCompound collectionsEnergy minimization with AMMOS and the Tff andMMFF94s force fields was carried out initially on 4 smallmolecules, taken from several X-ray protein-ligand com-plexes in the Protein Data Bank (PDB) [54], namely:raloxifene(an inhibitor of estrogen receptor, PDB code4-(n-acetylamino)-3-[n-(2-ethylbu-1err);tanoylamino)]benzoic acid (FDI) (an inhibitor of neu-raminidase, PDB code 1b9s);thymidine(an inhibitor ofthymidine kinase, PDB code 1kim); thieno [3,2-b]pyrid-ine-2-sulfonic acid [2-oxo-1-(1h-pyrrolo[2,3-c]pyridin-2-ylmethyl)- pyrrolidin-3-yl]-amide (PR2) (an inhibitor ofcoagulation factor X, PDB code 1f0r).To test the performance of AMMOS on a large number ofsmall organic molecules we used the ChemBridge diver-sity set [55]. Our decoy library contained 37970 mole-cules after ADME/Tox filtering with the program FilterPage 4 of 15(page number not for citation purposes)BMC Bioinformatics2008,9:438http://www.biomedcentral.com/1471-2105/9/4381.0.2 [43]. We merged 20 known active compounds forER and NA (with activities ranging from micromolar tonanomolar) and a number of rotatable bonds rangingfrom 4 to15) to the decoy collection. Some of the activecompounds were taken from the PDB protein-ligandstructures: 2 for ER (1err, 3ert; resolution 2.60 Å, 1.90 Å,respectively) and 10 for NA (1inf, 1inv, 1ivb, 1vcj, 1b9s,1b9t, 1b9v, 1a4g, 1f8b, 2qwk; resolution 2.35–2.50 Å).When a ligand could not be extracted from the PDB, it wasrebuilt from the literature [56]. All these active inhibitorswere added to the decoy library, all in SMILES format. Theresulting chemical library was transformed in single 3Dconformer and saved in mol2 format using the programOMEGA 2.0 [43]. The multiconformer states were thengenerated by our program Multiconf-DOCK [47] applyingan energy window of 25 kcal.mol-1and a diversity thresh-old of 1 Å RMSD. A maximum of 50 conformers were gen-erated for each molecule. These values represent anappropriate balance between speed and accuracy accord-ing to the recent studies [47,57].Protein targetsThe performance of AMMOS for post-processing of thepre-docked protein-ligand complexes was validated ontwo protein targets. We selected ER with a closed andhydrophobic pocket and NA with an open and highlypolar binding site. These two proteins present a bindingsite with very different degrees of burial (75.4% for ER and30.5% for NA) and polarity (25% for ER and 65% for NA)(see for details [20]). We took the co-crystallized struc-tures with best resolution among all retrieved protein-lig-and complexes (PDB code 3ert, resolution 1.90 Å for ERand PDB code 1b9v, resolution 2.35 Å for NA). All boundwater molecules and crystallized ligands were removedfrom the binding sites. Hydrogens were assigned using theprogram InsightII [46].and complexes: from fully flexible minimization of thewhole protein-ligand complex (case 1) to a flexible ligandin a rigid receptor (case 5) (see Fig. 1). In all situations, theligand atoms are free to move. Our minimization stepwith AMMP applied on the protein-ligand complexes forthe selected protein flexibility case involves 2×500 itera-tions with conjugate gradient optimization. The advancesuser can select any minimization method available inAMMP and specify the minimization parameters (i.e.number of iterations, convergence etc.). Finally, all theminimized conformers are scored by the AMMP mini-mized interaction protein-ligand energy [28] and re-ranked.Implementation of AMMOSThe package AMMOS consists of several programs devel-oped in C and Python, and makes use of the open sourceprograms AMMP and PREAMMP. One AMMOS routine(written in C) ensures the transformation of the input files(PDB for the protein and mol2 for the ligands) to a spe-cific ammp format required by AMMP and at the end ofthe process generate reversely, the protein in PDB formatand the small molecules in mol2 format. The automatiza-tion of the procedure described in Figure 1 for a largenumber of ligands in a single or multiple conformer stateis accomplished via a Python script. Five different cases(scripts written in C) have been elaborated for the selec-tion of the active/inactive atoms in the protein, while, inall situations, the ligand atoms are flexible:case 1:allatoms of the protein are active (a fully flexible minimiza-tion);case 2:all atoms of the protein side chains are active;case 3:all atoms of the protein inside a sphere around theligand are active;case 4:all atoms of the protein sidechains inside a sphere around the ligand are active;case 5:the whole protein is rigid. After processing of the wholepre-docked compound collection, the following resultsare saved in a subdirectory namedOUTPUT:(i) the coor-dinates of all conformers after minimization, (ii) the coor-dinates of the flexible part of the protein afterminimization, (iii) a file with warnings (if any), and (iv)the interaction energy between the protein and the ligands(external energy), the internal energies of the ligands orthe protein, and total energy (including internal andexternal terms) before and after minimization. Finally are-ranking step takes place based on the computed exter-nal energies of all minimized ligands. The method selectsthe best conformer among the multiple ligand conform-ers. The complete automatic procedure could be run foreither sp4 or sp5 force field. Because sp5 was not explicitlyavailable in AMMP, we created it using the extended andparameterized list of atom types and electronegativity val-ues available in reference [39].ResultsAlgorithm of AMMOSAMMOS drives a fully automatic procedure for minimiza-tion of protein-ligand complexes in a situation where thecompounds are pre-docked in the binding site by anydocking engine. AMMOS parameters are optimized suchas to handle relatively large docked compound collec-tions. Figure 1 illustrates all the different steps andrequired inputs, preparation, minimization of protein-lig-and complexes with different degrees of protein flexibil-ity, as well as final ranking of the ligands according to theminimized interaction energy between the ligands andthe receptor. Firstly, AMMOS employs the programPREAMMP to convert the input protein and ligand files toAMMP format. Next, AMMP autolink is run to search forincomplete amino-acid residues and to finally link all theresidues after corrections. AMMOS allows users to selectone among five different solutions to handle protein-lig-Another application of AMMOS is for the minimizationof a large databank of chemical compounds in thePage 5 of 15(page number not for citation purposes)

zanotowane.pl

doc.pisz.pl

pdf.pisz.pl

hannaeva.xlx.pl