Computational protocol: Novel phytochemical–antibiotic conjugates as multitarget inhibitors of Pseudomononas aeruginosa GyrB/ParE and DHFR

Similar protocols

Protocol publication

[…] Marvin sketch (ChemAxon, Budapest, Hungary) was used to create the SMILES (.smi) format for the designed inhibitors. Open Babel online interface (http://openbabel.sourceforge.net/), through shell script and using .smi as input, was used to generate primary 3D structures, which were then saved in SDF file format. The structures were subjected to conformational search calculations in gas phase using the Merck molecular force field method in the SPARTAN’10 program (WaveFunction, Inc, Irvine, CA, USA). The minimum energy conformation was chosen for geometry-optimization using the Austin model 1 (AM1) semiempirical algorithm. The optimized structures were reoptimized using the density functional theory hybrid, B3LYP/6-31G** basis set functional. The stereoelectronic properties such as HOMO (highest occupied molecular orbital) and LUMO (lowest unoccupied molecular orbital) energies, band gap (HOMO–LUMO), MEP, dipole moment and polarizability, physicochemical properties (molecular volume, polar surface area [PSA]), and absolute energies were calculated from single-point energy calculations on the complete geometry-optimized structures of the designed compounds using B3LYP/6-31G** basis set. In addition, percentage of absorption (%ABS) was calculated by using the following equation,based on the computed PSA values.To obtain additional insight on the toxicity risks (mutagenicity, tumorigenicity, irritation, and reproductive effects) and physicochemical properties (molecular weight, cLogP, solubility, PSA, number of hydrogen bond donors and acceptors, number of rotatable bonds, drug-likeness, and overall drug score) of the designed hybrid compounds, GyrB/ParE inhibitors (cyclothialidine and novobiocin) and DHFR inhibitors (methotrexate and trimethoprim), phytochemicals (protocatechuic acid and gallic acid), and sulfonamides (sulfamethoxazole and sulfadiazine), were calculated using OSIRIS property explorer and Molinspiration property prediction toolkit (Molinspiration Cheminformatics, Slovensky Grob, Slovak Republic) (http://www.molinspiration.com/cgi-bin/properties). The toxicity risk assessment indicated potential risks in the drawn structure concerning the risk category specified, and the process relied on the precomputed set of structural fragments that give toxicity alerts if they were encountered in the drawn 2D structures. Properties with high risks of undesired effects, like mutagenicity or a poor intestinal absorption, were shown as a negative sign, whereas a positive sign indicated drug-conform behavior. [...] The MolDock docking engine in Molegro virtual docker (CLC bio, Aarhus N, Denmark) was used for studying the binding modes of the designed hybrid compounds in the active site cavities of the modeled P. aeruginosa GyrB/ParE and DHFR enzymes (see Supplementary materials), which was based on a new heuristic search algorithm (MolDock score) that combines differential evolution with a cavity prediction algorithm. In our docking experiments, a MolDock grid scoring function using template docking with default values: −500 overall strength and 0.4 Å energy grid resolution was used to evaluate the energy between the ligand and the target enzyme. Grid resolution, number of runs, population size, maximum iterations, pose generation energy threshold, simplex evolution max steps, and neighbor distance factor were set as 0.30 Å, 20, 50, 1500, 100, 1.00 for each run, respectively, using the MolDock SE algorithm. The ligands from the crystal structures of Escherichia coli GyrB/ParE and DHFR were transferred into the workspace, keeping the orientation as a control and were kept as the reference ligand. The complete geometry-optimized structures of the hybrid compounds and the generated protein homologs were also transferred, and hydrogen molecules were added to both ligands and protein molecules using the preparation wizard in the Molegro workspace. During import of the 3D structures of the ligands, charges and bond orders were assigned, the torsional angle of the 3D structures was also determined, and all acyclic single bonds were set as flexible. Binding sites in the electrostatic surface of the protein were identified using the grid-based cavity prediction algorithm. A total of five cavities were detected, the prepositioned reference ligand in the active site cavity was identified, and the docking was constrained to the predicted active site cavity. Multiple poses were returned for each run with the root mean square deviation (RMSD) threshold set to 1.00 Å. The pose with the highest rerank and MolDock score was retained in the workspace for detailed evaluation of the ligand binding at the active site cavity. The rerank score uses a weighted combination of the terms used by the MolDock score mixed with a few additional terms (the rerank score includes the steric terms which are Lennard–Jones approximations to the steric energy; the MolDock score uses a piece-wise linear potential to approximate the steric energy). The rerank scoring function improved the docking accuracy by identifying the most promising docking solution from the solutions obtained by the MolDock docking algorithm. The rerank score provided an estimate of the strength of the interaction. It was not calibrated in chemical units, and it did not take complex contributions such as entropy into account. Even though the rerank score might be successful in ranking different poses of the same ligand, it might be less successful in ranking poses of different ligands. Along with both MolDock and reranking scores, we also predicted binding affinities using a calibrated model that is included in the Molegro virtual docker. The binding affinity model was trained using a data set of more than 200 structurally diverse complexes from Protein Data Bank (PDB) with known binding affinities. Hence, in our docking experiments we used this recommended strategy of ranking the docking results by their rerank scores and subsequently the binding affinity measure to get high ranked poses. The validation of the docking protocol was carried out by redocking the imported reference ligands from their respective experimental PDB structures in the predicted active site cavity of the model using the RMSD measure. The RMSD (Å) values for the redocked inhibitor methotrexate (MTX) and substrate NADP in the modeled P. aeruginosa DHFR active site were 0.3968 Å and 0.5794 Å, respectively, wherein the substrate, ADPNP (5′–adenylyl-β-γ-imidodiphosphate, the nonhydrolyzable analog of ATP) in the P. aeruginosa GyrB/ParE were 0.2529 Å/0.2052 Å, indicating high similarity between predicted and experimental binding mode. Furthermore, the binding conformations of ligands in the enzyme active site cavity were analyzed by visualizing the H-bond and electrostatic interactions formed between the ligands and the active site residues. [...] The molecular dynamics simulations for the corresponding substrate bound enzymes and top ranked docked inhibitor in both the enzymes were performed with the YASARA dynamics package version 10 (YASARA Biosciences, Vienna, Austria). A periodic simulation cell of at least 20 Å larger than the protein was employed with explicit solvent around the P. aeruginosa enzyme complexes (I–GyrB, IV–DHFR VII–ParE, II–ADPNP–GyrB, V–MTX-NADP-DHFR, VIII–ADPNP-ParE), and top ranked hybrid compound-bound enzyme complexes (III-phyto-drug [PD]_2a-GyrB, VI–PD_2a-NADP-DHFR, and IX–PD_2a-ParE) were energy minimized to correct the covalent geometry and to remove bumps. The AMBER03 force field was used with Van der Waals pairs cut-off distance 7.86 Å and long range statistics algorithm calculated using the Particle Mesh Ewald method. The automatic force field parameter assignment was carried out by AutoSMILES, which assigned pH-dependent fractional bond orders and protonation patterns. This was followed by geometry-optimization of the structures with the COSMO solvation model using semiempirical AM1 Mulliken point charge calculations and assignment of AM1BCC (AM1 bond charge correction), atom, and bond types, with further refinement using known restrained electrostatic potential charges. And finally, GAFF (General AMBER force field) atom types and remaining force field parameters were assigned. The hydrogen bonding network of the complexes were optimized by a method that positions polar hydrogen atoms in protein structures by optimizing the total hydrogen bond energy. The simulation cell was filled with TIP3 water, and Na+ and Cl− counter ions were added at the most favorable positions to neutralize the ion strength in the cell. A few water molecules were deleted to readjust the solvent density to 0.997 g/mL, and the pressure was controlled by rescaling the simulation box dimensions to maintain the water density. The pKa values were predicted using the empirical pKa prediction equation that is approximated as a function of electrostatic potential, hydrogen bonds, and accessible surface area by Ewald summation implemented in YASARA. A short molecular dynamics simulation was run on the solvent only, and the entire system was energy minimized by steepest descent minimization to remove conformational stress. This was followed by simulated annealing minimization using a time step of 2 femtoseconds and atom velocities scaled down by 0.9 every 10th step until convergence (energy improved by less than 0.1% during 200 steps). The molecular dynamics simulation was then initiated, with temperature fixed at 298 K and multiple time steps set for intramolecular and intermolecular forces at 1.33 femtoseconds and 4 femtoseconds, respectively. Each complex was subjected to 3000 picoseconds (ps) of molecular dynamics simulations, and snapshots were saved every 1 ps for data analysis. The trajectories of the molecular dynamics simulations were analyzed for the equilibrium stability by measuring the RMSD of the complexes and the root mean square fluctuation (RMSF) of residues around the ligand-binding active sites. […]

Pipeline specifications

Software tools Open Babel, Molinspiration, MVD, YASARA Dynamics, YASARA
Applications Drug design, Protein interaction analysis
Organisms Pseudomonas aeruginosa
Diseases Pseudomonas Infections
Chemicals Adenosine Triphosphate, Folic Acid, Hydrogen