1 #!/bin/env python
2 #
3 # File: VinaPerformDocking.py
4 # Author: Manish Sud <msud@san.rr.com>
5 #
6 # Acknowledgments: Diogo Santos-Martins and Stefano Forli
7 #
8 # Copyright (C) 2024 Manish Sud. All rights reserved.
9 #
10 # The functionality available in this script is implemented using AutoDockVina
11 # and Meeko, open source software packages for docking, and RDKit, an open
12 # source toolkit for cheminformatics developed by Greg Landrum.
13 #
14 # This file is part of MayaChemTools.
15 #
16 # MayaChemTools is free software; you can redistribute it and/or modify it under
17 # the terms of the GNU Lesser General Public License as published by the Free
18 # Software Foundation; either version 3 of the License, or (at your option) any
19 # later version.
20 #
21 # MayaChemTools is distributed in the hope that it will be useful, but without
22 # any warranty; without even the implied warranty of merchantability of fitness
23 # for a particular purpose. See the GNU Lesser General Public License for more
24 # details.
25 #
26 # You should have received a copy of the GNU Lesser General Public License
27 # along with MayaChemTools; if not, see <http://www.gnu.org/licenses/> or
28 # write to the Free Software Foundation Inc., 59 Temple Place, Suite 330,
29 # Boston, MA, 02111-1307, USA.
30 #
31
32 from __future__ import print_function
33
34 # Add local python path to the global path and import standard library modules...
35 import os
36 import sys; sys.path.insert(0, os.path.join(os.path.dirname(sys.argv[0]), "..", "lib", "Python"))
37 import time
38 import re
39 import shutil
40 import tempfile
41 import json
42 import glob
43 import multiprocessing as mp
44
45 # AutoDock Vina imports...
46 try:
47 from vina import Vina
48 from vina import __version__ as vinaVersion
49 except ImportError as ErrMsg:
50 sys.stderr.write("\nFailed to import AutoDock Vina module/package: %s\n" % ErrMsg)
51 sys.stderr.write("Check/update your Vina environment and try again.\n\n")
52 sys.exit(1)
53
54 # AutoDock Meeko imports...
55 try:
56 from meeko import __version__ as meekoVersion
57 from meeko import MoleculePreparation
58 from meeko import PDBQTMolecule
59 from meeko import RDKitMolCreate
60 from meeko import PDBQTWriterLegacy
61 except ImportError as ErrMsg:
62 sys.stderr.write("\nFailed to import AutoDock Meeko module/package: %s\n" % ErrMsg)
63 sys.stderr.write("Check/update your Meeko environment and try again.\n\n")
64 sys.exit(1)
65
66 # RDKit imports...
67 try:
68 from rdkit import rdBase
69 from rdkit import Chem
70 from rdkit.Chem import AllChem
71 from rdkit.Chem import rdMolTransforms
72 except ImportError as ErrMsg:
73 sys.stderr.write("\nFailed to import RDKit module/package: %s\n" % ErrMsg)
74 sys.stderr.write("Check/update your RDKit environment and try again.\n\n")
75 sys.exit(1)
76
77 # MayaChemTools imports...
78 try:
79 from docopt import docopt
80 import MiscUtil
81 import RDKitUtil
82 except ImportError as ErrMsg:
83 sys.stderr.write("\nFailed to import MayaChemTools module/package: %s\n" % ErrMsg)
84 sys.stderr.write("Check/update your MayaChemTools environment and try again.\n\n")
85 sys.exit(1)
86
87 ScriptName = os.path.basename(sys.argv[0])
88 Options = {}
89 OptionsInfo = {}
90
91 def main():
92 """Start execution of the script."""
93
94 MiscUtil.PrintInfo("\n%s (Vina v%s; Meeko v%s; RDKit v%s; MayaChemTools v%s; %s): Starting...\n" % (ScriptName, vinaVersion, meekoVersion, rdBase.rdkitVersion, MiscUtil.GetMayaChemToolsVersion(), time.asctime()))
95
96 (WallClockTime, ProcessorTime) = MiscUtil.GetWallClockAndProcessorTime()
97
98 # Retrieve command line arguments and options...
99 RetrieveOptions()
100
101 # Process and validate command line arguments and options...
102 ProcessOptions()
103
104 # Perform actions required by the script...
105 PerformDocking()
106
107 MiscUtil.PrintInfo("\n%s: Done...\n" % ScriptName)
108 MiscUtil.PrintInfo("Total time: %s" % MiscUtil.GetFormattedElapsedTime(WallClockTime, ProcessorTime))
109
110 def PerformDocking():
111 """Perform docking."""
112
113 # Setup a molecule reader...
114 MiscUtil.PrintInfo("\nProcessing file %s..." % OptionsInfo["Infile"])
115 Mols = RDKitUtil.ReadMolecules(OptionsInfo["Infile"], **OptionsInfo["InfileParams"])
116
117 # Set up molecule writers...
118 Writer, WriterFlexRes = SetupWriters()
119 if WriterFlexRes is None:
120 MiscUtil.PrintInfo("Generating file %s..." % OptionsInfo["Outfile"])
121 else:
122 MiscUtil.PrintInfo("Generating files %s and %s..." % (OptionsInfo["Outfile"], OptionsInfo["OutfileFlexRes"]))
123
124 MolCount, ValidMolCount, DockingFailedCount = ProcessMolecules(Mols, Writer, WriterFlexRes)
125
126 CloseWriters(Writer, WriterFlexRes)
127
128 MiscUtil.PrintInfo("\nTotal number of molecules: %d" % MolCount)
129 MiscUtil.PrintInfo("Number of valid molecules: %d" % ValidMolCount)
130 MiscUtil.PrintInfo("Number of molecules failed during docking: %d" % DockingFailedCount)
131 MiscUtil.PrintInfo("Number of ignored molecules: %d" % (MolCount - ValidMolCount + DockingFailedCount))
132
133 def ProcessMolecules(Mols, Writer, WriterFlexRes):
134 """Process and dock molecules."""
135
136 if OptionsInfo["MPMode"]:
137 return ProcessMoleculesUsingMultipleProcesses(Mols, Writer, WriterFlexRes)
138 else:
139 return ProcessMoleculesUsingSingleProcess(Mols, Writer, WriterFlexRes)
140
141 def ProcessMoleculesUsingSingleProcess(Mols, Writer, WriterFlexRes):
142 """Process and dock molecules using a single process."""
143
144 VinaHandle = InitializeVina(OptionsInfo["QuietMode"])
145 MolPrepHandle = InitializeMeekoMolPrepration(OptionsInfo["QuietMode"])
146 if not OptionsInfo["QuietMode"]:
147 MiscUtil.PrintInfo("\nVina configuration:\n%s" % VinaHandle)
148
149 MiscUtil.PrintInfo("\nDocking molecules (Mode: %s)..." % OptionsInfo["Mode"])
150
151 (MolCount, ValidMolCount, DockingFailedCount) = [0] * 3
152 for Mol in Mols:
153 MolCount += 1
154
155 if Mol is None:
156 continue
157
158 if not CheckAndValidateMolecule(Mol, MolCount):
159 continue
160
161 ValidMolCount += 1
162 CalcStatus, PoseMol, PoseMolEnergies, PoseMolFlexRes = DockMolecule(VinaHandle, MolPrepHandle, Mol, MolCount)
163
164 if not CalcStatus:
165 if not OptionsInfo["QuietMode"]:
166 MiscUtil.PrintWarning("Failed to dock molecule %s" % RDKitUtil.GetMolName(Mol, MolCount))
167
168 DockingFailedCount += 1
169 continue
170
171 WriteMolPoses(Writer, Mol, MolCount, PoseMol, PoseMolEnergies, WriterFlexRes, PoseMolFlexRes)
172
173 return (MolCount, ValidMolCount, DockingFailedCount)
174
175 def ProcessMoleculesUsingMultipleProcesses(Mols, Writer, WriterFlexRes):
176 """Process and calculate energy of molecules using multiprocessing."""
177
178 MiscUtil.PrintInfo("\nDocking molecules using multiprocessing...")
179 MiscUtil.PrintInfo("\nDocking molecules (Mode: %s)..." % OptionsInfo["Mode"])
180
181 MPParams = OptionsInfo["MPParams"]
182
183 # Setup data for initializing a worker process...
184 InitializeWorkerProcessArgs = (MiscUtil.ObjectToBase64EncodedString(Options), MiscUtil.ObjectToBase64EncodedString(OptionsInfo))
185
186 # Setup a encoded mols data iterable for a worker process...
187 WorkerProcessDataIterable = RDKitUtil.GenerateBase64EncodedMolStrings(Mols)
188
189 # Setup process pool along with data initialization for each process...
190 if not OptionsInfo["QuietMode"]:
191 MiscUtil.PrintInfo("\nConfiguring multiprocessing using %s method..." % ("mp.Pool.imap()" if re.match("^Lazy$", MPParams["InputDataMode"], re.I) else "mp.Pool.map()"))
192 MiscUtil.PrintInfo("NumProcesses: %s; InputDataMode: %s; ChunkSize: %s\n" % (MPParams["NumProcesses"], MPParams["InputDataMode"], ("automatic" if MPParams["ChunkSize"] is None else MPParams["ChunkSize"])))
193
194 ProcessPool = mp.Pool(MPParams["NumProcesses"], InitializeWorkerProcess, InitializeWorkerProcessArgs)
195
196 # Start processing...
197 if re.match("^Lazy$", MPParams["InputDataMode"], re.I):
198 Results = ProcessPool.imap(WorkerProcess, WorkerProcessDataIterable, MPParams["ChunkSize"])
199 elif re.match("^InMemory$", MPParams["InputDataMode"], re.I):
200 Results = ProcessPool.map(WorkerProcess, WorkerProcessDataIterable, MPParams["ChunkSize"])
201 else:
202 MiscUtil.PrintError("The value, %s, specified for \"--inputDataMode\" is not supported." % (MPParams["InputDataMode"]))
203
204 (MolCount, ValidMolCount, DockingFailedCount) = [0] * 3
205 for Result in Results:
206 MolCount += 1
207
208 MolIndex, EncodedMol, CalcStatus, EncodedPoseMol, PoseMolEnergies, EncodedPoseMolFlexRes = Result
209
210 if EncodedMol is None:
211 continue
212
213 ValidMolCount += 1
214
215 if not CalcStatus:
216 if not OptionsInfo["QuietMode"]:
217 MolName = RDKitUtil.GetMolName(Mol, MolCount)
218 MiscUtil.PrintWarning("Failed to dock molecule %s" % MolName)
219
220 DockingFailedCount += 1
221 continue
222
223 Mol = RDKitUtil.MolFromBase64EncodedMolString(EncodedMol)
224 PoseMol = None if EncodedPoseMol is None else RDKitUtil.MolFromBase64EncodedMolString(EncodedPoseMol)
225 PoseMolFlexRes = None if EncodedPoseMolFlexRes is None else RDKitUtil.MolFromBase64EncodedMolString(EncodedPoseMolFlexRes)
226
227 WriteMolPoses(Writer, Mol, MolCount, PoseMol, PoseMolEnergies, WriterFlexRes, PoseMolFlexRes)
228
229 return (MolCount, ValidMolCount, DockingFailedCount)
230
231 def InitializeWorkerProcess(*EncodedArgs):
232 """Initialize data for a worker process."""
233
234 global Options, OptionsInfo
235
236 if not OptionsInfo["QuietMode"]:
237 MiscUtil.PrintInfo("Starting process (PID: %s)..." % os.getpid())
238
239 # Decode Options and OptionInfo...
240 Options = MiscUtil.ObjectFromBase64EncodedString(EncodedArgs[0])
241 OptionsInfo = MiscUtil.ObjectFromBase64EncodedString(EncodedArgs[1])
242
243 # Initialize in a quiet mode...
244 OptionsInfo["VinaHandle"] = InitializeVina(True)
245 OptionsInfo["MolPrepHandle"] = InitializeMeekoMolPrepration(True)
246
247 def WorkerProcess(EncodedMolInfo):
248 """Process data for a worker process."""
249
250 MolIndex, EncodedMol = EncodedMolInfo
251
252 CalcStatus, PoseMol, PoseMolEnergies, PoseMolFlexRes = [False, None, None, None]
253
254 if EncodedMol is None:
255 return [MolIndex, None, CalcStatus, PoseMol, Energies, PoseMolFlexRes]
256
257 Mol = RDKitUtil.MolFromBase64EncodedMolString(EncodedMol)
258 MolCount = MolIndex + 1
259
260 if not CheckAndValidateMolecule(Mol, MolCount):
261 return [MolIndex, None, CalcStatus, PoseMol, Energies, PoseMolFlexRes]
262
263 CalcStatus, PoseMol, PoseMolEnergies, PoseMolFlexRes = DockMolecule(OptionsInfo["VinaHandle"], OptionsInfo["MolPrepHandle"], Mol, MolCount)
264
265 EncodedPoseMol = None if PoseMol is None else RDKitUtil.MolToBase64EncodedMolString(PoseMol, PropertyPickleFlags = Chem.PropertyPickleOptions.MolProps | Chem.PropertyPickleOptions.PrivateProps)
266
267 EncodedPoseMolFlexRes = None if PoseMolFlexRes is None else RDKitUtil.MolToBase64EncodedMolString(PoseMolFlexRes, PropertyPickleFlags = Chem.PropertyPickleOptions.MolProps | Chem.PropertyPickleOptions.PrivateProps)
268
269 return [MolIndex, EncodedMol, CalcStatus, EncodedPoseMol, PoseMolEnergies, EncodedPoseMolFlexRes]
270
271 def DockMolecule(VinaHandle, MolPrepHandle, Mol, MolNum = None):
272 """Dock molecule."""
273
274 Status, PoseMol, PoseMolEnergies, PoseMolFlexRes = [False, None, None, None]
275
276 if OptionsInfo["VinaVerbosity"] != 0:
277 MiscUtil.PrintInfo("\nProcessing molecule %s..." % RDKitUtil.GetMolName(Mol, MolNum))
278
279 PDBQTMolStr = PrepareMolecule(MolPrepHandle, Mol)
280 if PDBQTMolStr is None:
281 return (Status, PoseMol, PoseMolEnergies, PoseMolFlexRes)
282 VinaHandle.set_ligand_from_string(PDBQTMolStr)
283
284 if OptionsInfo["ScoreOnlyMode"]:
285 return ProcessMoleculeForScoreOnlyMode(VinaHandle, Mol)
286 elif OptionsInfo["LocalOptimizationOnlyMode"]:
287 return ProcessMoleculeForLocalOptimizationOnlyMode(VinaHandle, Mol)
288 elif OptionsInfo["DockMode"]:
289 return ProcessMoleculeForDockMode(VinaHandle, Mol)
290 else:
291 MiscUtil.PrintError("The value specified, %s, for option \"-m, --mode\" is not valid. Supported values: Dock LocalOptimizationOnly ScoreOnly" % OptionsInfo["Mode"])
292
293 return (Status, PoseMol, PoseMolEnergies, PoseMolFlexRes)
294
295 def ProcessMoleculeForScoreOnlyMode(VinaHandle, Mol):
296 """Score molecule without performing any docking."""
297
298 Status, PoseMol, Energies, PoseMolFlexRes = [False, None, None, None]
299
300 # Score molecule...
301 try:
302 Energies = VinaHandle.score()
303 except Exception as ErrMsg:
304 if not OptionsInfo["QuietMode"]:
305 MiscUtil.PrintWarning("Failed to score molecule:\n%s\n" % ErrMsg)
306 return (False, None, None, None)
307
308 if len(Energies) == 0:
309 return (False, None, None, None)
310
311 Status = True
312 PoseMol = Mol
313
314 return (Status, PoseMol, Energies, PoseMolFlexRes)
315
316 def ProcessMoleculeForLocalOptimizationOnlyMode(VinaHandle, Mol):
317 """Score molecule after a local optimization wthout any docking."""
318
319 Status, PoseMol, Energies, PoseMolFlexRes = [False, None, None, None]
320
321 # Optimize and score molecule...
322 try:
323 Energies = VinaHandle.optimize()
324 except Exception as ErrMsg:
325 if not OptionsInfo["QuietMode"]:
326 MiscUtil.PrintWarning("Failed to perform local optimization:\n%s\n" % ErrMsg)
327 return (False, None, None, None)
328
329 if len(Energies) == 0:
330 return (False, None, None, None)
331
332 # Write optimize pose to a temporray file and retrieve the PDBQT string for the pose...
333 (_, TmpFile) = tempfile.mkstemp(suffix = ".pdbqt", prefix = "VinaOptimize_", text = True)
334 VinaHandle.write_pose(TmpFile, overwrite = True)
335
336 TmpFH = open(TmpFile, "r")
337 PosePDBQTOutputStr = TmpFH.read()
338 TmpFH.close()
339
340 os.remove(TmpFile)
341
342 # Setup a mol containing optimized pose...
343 try:
344 PosePDBQTOutputMol = PDBQTMolecule(PosePDBQTOutputStr)
345 PoseMols = RDKitMolCreate.from_pdbqt_mol(PosePDBQTOutputMol)
346 except Exception as ErrMsg:
347 if not OptionsInfo["QuietMode"]:
348 MiscUtil.PrintWarning("Failed to retrieve optimize pose:\n%s\n" % ErrMsg)
349 return (False, None, None, None)
350
351 if len(PoseMols) == 0:
352 return (False, None, None, None)
353
354 Status = True
355 PoseMol = PoseMols[0]
356
357 return (Status, PoseMol, Energies, PoseMolFlexRes)
358
359 def ProcessMoleculeForDockMode(VinaHandle, Mol):
360 """Dock and score molecule."""
361
362 Status, PoseMol, Energies, PoseMolFlexRes = [False, None, None, None]
363
364 # Dock molecule...
365 try:
366 VinaHandle.dock(exhaustiveness = OptionsInfo["Exhaustiveness"] , n_poses = OptionsInfo["NumPoses"], min_rmsd = OptionsInfo["MinRMSD"], max_evals = OptionsInfo["MaxEvaluations"])
367 except Exception as ErrMsg:
368 if not OptionsInfo["QuietMode"]:
369 MiscUtil.PrintWarning("Failed to dock molecule:\n%s\n" % ErrMsg)
370 return (False, None, None, None)
371
372 # Retrieve poses...
373 try:
374 PosePDBQTOutputStr = VinaHandle.poses(n_poses = OptionsInfo["NumPoses"], energy_range = OptionsInfo["EnergyRange"], coordinates_only = False)
375 except Exception as ErrMsg:
376 if not OptionsInfo["QuietMode"]:
377 MiscUtil.PrintWarning("Failed to retrieve docked poses:\n%s\n" % ErrMsg)
378 return (False, None, None, None)
379
380 PoseMol, PoseMolFlexRes = SetupPoseMols(PosePDBQTOutputStr)
381 if PoseMol is None:
382 return (False, None, None, None)
383
384 # Retrieve energies...
385 try:
386 Energies = VinaHandle.energies(n_poses = OptionsInfo["NumPoses"], energy_range = OptionsInfo["EnergyRange"])
387 except Exception as ErrMsg:
388 if not OptionsInfo["QuietMode"]:
389 MiscUtil.PrintWarning("Failed to retrieve energies for docked poses:\n%s\n" % ErrMsg)
390 return (False, None, None, None)
391
392 if len(Energies) == 0:
393 return (False, None, None, None)
394
395 Status = True
396
397 return (Status, PoseMol, Energies, PoseMolFlexRes)
398
399 def SetupPoseMols(PosePDBQTOutputStr):
400 """Process PDBQT Vina poses string to setup pose mols for the docked
401 molecule and any flexible residues."""
402
403 PoseMol, PoseMolFlexRes = [None, None]
404
405 # Setup a mol containing poses as conformers...
406 try:
407 PosePDBQTOutputMol = PDBQTMolecule(PosePDBQTOutputStr)
408 PoseMols = RDKitMolCreate.from_pdbqt_mol(PosePDBQTOutputMol)
409 except Exception as ErrMsg:
410 if not OptionsInfo["QuietMode"]:
411 MiscUtil.PrintWarning("Failed to retrieve docked poses:\n%s\n" % ErrMsg)
412 return (None, None)
413
414 if len(PoseMols) == 0:
415 return (None, None)
416
417 # First mol is the pose for docked molecule...
418 PoseMol = PoseMols.pop(0)
419
420 # Collect pose mols foe flexible side chain residues...
421 PoseMolsFlexRes = []
422 for Mol in PoseMols:
423 if not Mol.HasProp("meeko"):
424 continue
425 MeekoPropMap = json.loads(Mol.GetProp("meeko"))
426 if type(MeekoPropMap) is dict:
427 if "is_sidechain" in MeekoPropMap:
428 if MeekoPropMap["is_sidechain"]:
429 PoseMolsFlexRes.append(Mol)
430
431 # Combine pose mols for flexible side chain residues into a single pose mol...
432 if len(PoseMolsFlexRes):
433 for Mol in PoseMolsFlexRes:
434 if PoseMolFlexRes is None:
435 PoseMolFlexRes = Mol
436 else:
437 PoseMolFlexRes = Chem.CombineMols(PoseMolFlexRes, Mol)
438
439 if PoseMolFlexRes is not None:
440 PoseMolConfCount = PoseMol.GetNumConformers()
441 PoseMolFlexResConfCount = PoseMolFlexRes.GetNumConformers()
442 if PoseMolConfCount != PoseMolFlexResConfCount:
443 if not OptionsInfo["QuietMode"]:
444 MiscUtil.PrintWarning("The number of poses, %s, for flexible residues doesn't match number of poses, %s, for docked molecule...\n" % (PoseMolFlexResConfCount, PoseMolConfCount))
445
446 return (PoseMol, PoseMolFlexRes)
447
448 def PrepareMolecule(MolPrepHandle, Mol):
449 """Prepare molecule for docking. """
450
451 try:
452 PreppedMols = MolPrepHandle.prepare(Mol)
453 except Exception as ErrMsg:
454 if not OptionsInfo["QuietMode"]:
455 MiscUtil.PrintWarning("Failed to prepare molecule for docking:\n%s\n" % ErrMsg)
456 return None
457
458 if len(PreppedMols) == 0:
459 return None
460
461 PreppedMol = PreppedMols[0]
462
463 # Setup PDBQT mole string...
464 try:
465 PDBQTMolStr, Status, ErrMsg = PDBQTWriterLegacy.write_string(PreppedMol)
466 except Exception as ExceptionErrMsg:
467 if not OptionsInfo["QuietMode"]:
468 MiscUtil.PrintWarning("Failed to prepare molecule for docking:\n%s\n" % ExceptionErrMsg)
469 return None
470
471 if not Status:
472 if not OptionsInfo["QuietMode"]:
473 MiscUtil.PrintWarning("Failed to prepare molecule for docking:\n%s\n" % ErrMsg)
474 return None
475
476 if MiscUtil.IsEmpty(PDBQTMolStr):
477 return None
478
479 return PDBQTMolStr
480
481 def InitializeVina(Quiet = False):
482 """Initialize AutoDock Vina. """
483
484 if not Quiet:
485 MiscUtil.PrintInfo("\nInitializing Vina...")
486
487 VinaHandle = Vina(sf_name = OptionsInfo["Forcefield"], cpu = OptionsInfo["NumThreads"], seed = OptionsInfo["RandomSeed"], no_refine = OptionsInfo["SkipRefinement"], verbosity = OptionsInfo["VinaVerbosity"])
488
489 SetupReceptor(VinaHandle, Quiet)
490 SetupForcefieldWeights(VinaHandle, Quiet)
491 SetupReceptorMaps(VinaHandle, Quiet)
492
493 return VinaHandle
494
495 def SetupReceptor(VinaHandle, Quiet = False):
496 """Setup receptor """
497
498 if OptionsInfo["UseReceptorFile"] and OptionsInfo["UseReceptorFlexFile"]:
499 VinaHandle.set_receptor(rigid_pdbqt_filename = OptionsInfo["ReceptorFile"], flex_pdbqt_filename = OptionsInfo["ReceptorFlexFile"])
500 elif OptionsInfo["UseReceptorFile"]:
501 VinaHandle.set_receptor(rigid_pdbqt_filename = OptionsInfo["ReceptorFile"], flex_pdbqt_filename = None)
502 elif OptionsInfo["UseReceptorFlexFile"]:
503 VinaHandle.set_receptor(rigid_pdbqt_filename = None, flex_pdbqt_filename = OptionsInfo["ReceptorFlexFile"])
504
505 def SetupForcefieldWeights(VinaHandle, Quiet = False):
506 """Setup forcefield weights. """
507
508 Weights = OptionsInfo["ForcefieldWeightParams"]
509 Forcefield = OptionsInfo["Forcefield"]
510 if not OptionsInfo["ForcefieldWeightParamsSpecified"]:
511 if not Quiet:
512 MiscUtil.PrintInfo("\nUsing default forcefield weights for %s..." % Forcefield)
513 return
514
515 if not Quiet:
516 MiscUtil.PrintInfo("\nSetting specified forcefield weights for %s..." % Forcefield)
517
518 if OptionsInfo["UseAD4Forcefield"]:
519 VinaHandle.set_weights([Weights["AD4Vdw"], Weights["AD4HydrogenBond"], Weights["AD4Electrostatic"], Weights["AD4Desolvation"], Weights["AD4GlueLinearAttraction"], Weights["AD4Rot"]])
520 elif OptionsInfo["UseVinaForcefield"]:
521 VinaHandle.set_weights([Weights["VinaGaussian1"], Weights["VinaGaussian2"], Weights["VinaRepulsion"], Weights["VinaHydrophobic"], Weights["VinaHydrogenBond"], Weights["VinaGlueLinearAttraction"], Weights["VinaRot"]])
522 elif OptionsInfo["UseVinardoForcefield"]:
523 VinaHandle.set_weights([Weights["VinardoGaussian1"], Weights["VinardoRepulsion"], Weights["VinardoHydrophobic"], Weights["VinardoHydrogenBond"], Weights["VinardoGlueLinearAttraction"], Weights["VinardoRot"]])
524 else:
525 MiscUtil.PrintError("The value specified, %s, for option \"-f, --forcefield\" is not valid. Supported values: AD4 Vina Vinardo" % Forcefield)
526
527 def SetupReceptorMaps(VinaHandle, Quiet = False):
528 """Setup receptor maps."""
529
530 if not Quiet:
531 MiscUtil.PrintInfo("\nSetting up receptor and maps for %s forcefield..." % OptionsInfo["Forcefield"])
532
533 if OptionsInfo["UseAD4Forcefield"]:
534 # Load maps for rigid portion of the receptor...
535 VinaHandle.load_maps(OptionsInfo["ReceptorMapsPrefix"])
536 elif OptionsInfo["UseVinaForcefield"] or OptionsInfo["UseVinardoForcefield"]:
537 # Load or compute maps for rigid portion of the receptor...
538 if OptionsInfo["UseReceptorMapsPrefix"]:
539 VinaHandle.load_maps(OptionsInfo["ReceptorMapsPrefix"])
540 else:
541 VinaHandle.compute_vina_maps(center = OptionsInfo["GridCenterList"], box_size = OptionsInfo["GridSizeList"], spacing = OptionsInfo["GridSpacing"], force_even_voxels = False)
542 else:
543 MiscUtil.PrintError("The value specified, %s, for option \"-f, --forcefield\" is not valid. Supported values: AD4 Vina Vinardo" % OptionsInfo["Forcefield"])
544
545 def InitializeMeekoMolPrepration(Quiet = False):
546 """Initialize meeko molecule prepration."""
547
548 if OptionsInfo["MergeHydrogens"]:
549 MolPrep = MoleculePreparation(merge_these_atom_types=("H",))
550 else:
551 MolPrep = MoleculePreparation(merge_these_atom_types="")
552
553 return MolPrep
554
555 def CheckAndValidateMolecule(Mol, MolCount = None):
556 """Validate molecule for docking."""
557
558 MolName = RDKitUtil.GetMolName(Mol, MolCount)
559 if RDKitUtil.IsMolEmpty(Mol):
560 if not OptionsInfo["QuietMode"]:
561 MiscUtil.PrintWarning("Ignoring empty molecule: %s\n" % MolName)
562 return False
563
564 if not OptionsInfo["ValidateMolecules"]:
565 return True
566
567 # Check for 3D flag...
568 if not Mol.GetConformer().Is3D():
569 if not OptionsInfo["QuietMode"]:
570 MiscUtil.PrintWarning("3D tag is not set for molecule: %s\n" % MolName)
571
572 # Check for missing hydrogens...
573 if RDKitUtil.AreHydrogensMissingInMolecule(Mol):
574 if not OptionsInfo["QuietMode"]:
575 MiscUtil.PrintWarning("Missing hydrogens in molecule: %s\n" % MolName)
576
577 return True
578
579 def WriteMolPoses(Writer, Mol, MolNum, PoseMol, Energies, WriterFlexRes = None, PoseMolFlexRes = None):
580 """Write molecule."""
581
582 if OptionsInfo["ScoreOnlyMode"]:
583 return WriteMolPoseForScoreOnlyMode(Writer, Mol, MolNum, PoseMol, Energies)
584 elif OptionsInfo["LocalOptimizationOnlyMode"]:
585 return WriteMolPoseForLocalOptimizationOnlyMode(Writer, Mol, MolNum, PoseMol, Energies)
586 elif OptionsInfo["DockMode"]:
587 return WriteMolPosesForDockMode(Writer, Mol, MolNum, PoseMol, Energies, WriterFlexRes, PoseMolFlexRes)
588 else:
589 MiscUtil.PrintError("The value specified, %s, for option \"-m, --mode\" is not valid. Supported values: Dock LocalOptimizationOnly ScoreOnly" % OptionsInfo["Mode"])
590
591 def WriteMolPoseForScoreOnlyMode(Writer, Mol, MolNum, PoseMol, Energies):
592 """Write out molecule and associated information for score only mode."""
593
594 ClearMeekoMolProperties(PoseMol)
595
596 MolName = RDKitUtil.GetMolName(Mol, MolNum)
597 PoseMol.SetProp("_Name", MolName)
598
599 SetupEnergyProperties(PoseMol, Energies)
600 Writer.write(PoseMol)
601
602 return
603
604 def WriteMolPoseForLocalOptimizationOnlyMode(Writer, Mol, MolNum, PoseMol, Energies):
605 """Write out molecule and associated information for score only mode."""
606
607 ClearMeekoMolProperties(PoseMol)
608
609 MolName = RDKitUtil.GetMolName(Mol, MolNum)
610 PoseMol.SetProp("_Name", MolName)
611
612 SetupEnergyProperties(PoseMol, Energies)
613 Writer.write(PoseMol)
614
615 def WriteMolPosesForDockMode(Writer, Mol, MolNum, PoseMol, Energies, WriterFlexRes = None, PoseMolFlexRes = None):
616 """Write out molecule and associated information for dock mode."""
617
618 MolName = RDKitUtil.GetMolName(Mol, MolNum)
619
620 # Write out poses for docked molecule...
621 ClearMeekoMolProperties(PoseMol)
622 for PoseMolConfIndex, PoseMolConf in enumerate(PoseMol.GetConformers()):
623 PoseMolName = "%s_Pose%s" % (MolName, (PoseMolConfIndex + 1))
624 PoseMol.SetProp("_Name", PoseMolName)
625
626 SetupEnergyProperties(PoseMol, Energies[PoseMolConfIndex])
627
628 Writer.write(PoseMol, confId = PoseMolConf.GetId())
629
630 # Write out poses for flexible reside side chains...
631 if WriterFlexRes is None or PoseMolFlexRes is None:
632 return
633
634 ClearMeekoMolProperties(PoseMolFlexRes)
635 for PoseMolFlexResConfIndex, PoseMolFlexResConf in enumerate(PoseMolFlexRes.GetConformers()):
636 PoseMolFlexResName = "%s_Flex_Receptor_Pose%s" % (MolName, (PoseMolFlexResConfIndex + 1))
637 PoseMolFlexRes.SetProp("_Name", PoseMolFlexResName)
638
639 WriterFlexRes.write(PoseMolFlexRes, confId = PoseMolFlexResConf.GetId())
640
641 def ClearMeekoMolProperties(Mol):
642 """Clear Meeko molecule properties. """
643
644 for PropName in ["meeko"]:
645 if Mol.HasProp(PropName):
646 Mol.ClearProp(PropName)
647
648 def SetupEnergyProperties(Mol, Energies):
649 """Setup energy properties. """
650
651 Precision = OptionsInfo["Precision"]
652
653 for Index, EnergyLabel in enumerate(OptionsInfo["EnergyLabelsList"]):
654 EnergyValueIndex = OptionsInfo["EnergyValueIndicesList"][Index]
655 EnergyValue = "%.*f" % (Precision, Energies[EnergyValueIndex])
656 Mol.SetProp(EnergyLabel, EnergyValue)
657
658 def SetupWriters():
659 """Setup writers for output files. """
660
661 Writer, WriterFlexRes = [None, None]
662
663 Writer = RDKitUtil.MoleculesWriter(OptionsInfo["Outfile"], **OptionsInfo["OutfileParams"])
664 if Writer is None:
665 MiscUtil.PrintError("Failed to setup a writer for output fie %s " % OptionsInfo["Outfile"])
666
667 if OptionsInfo["DockMode"] and OptionsInfo["UseReceptorFlexFile"]:
668 WriterFlexRes = RDKitUtil.MoleculesWriter(OptionsInfo["OutfileFlexRes"], **OptionsInfo["OutfileParams"])
669 if WriterFlexRes is None:
670 MiscUtil.PrintError("Failed to setup a writer for output fie %s " % OptionsInfo["OutfileFlexRes"])
671
672 return (Writer, WriterFlexRes)
673
674 def CloseWriters(Writer, WriterFlexRes):
675 """Close writers. """
676
677 if Writer is not None:
678 Writer.close()
679
680 if WriterFlexRes is not None:
681 WriterFlexRes.close()
682
683 def ComputeGridCenter(LigandFile):
684 """Compute grid center from ligand file."""
685
686 GridCenter = []
687
688 MiscUtil.PrintInfo("\nComputing grid center from reference ligand file %s..." % LigandFile)
689
690 Mols = RDKitUtil.ReadMolecules(LigandFile)
691 Mol = Mols[0]
692
693 Centroid = rdMolTransforms.ComputeCentroid(Mol.GetConformer())
694
695 GridCenter = [Centroid.x, Centroid.y, Centroid.z]
696
697 GridCenterFormatted = ["%.3f" % Value for Value in GridCenter]
698 MiscUtil.PrintInfo("GridCenter: %s" % (" ".join(GridCenterFormatted)))
699
700 return GridCenter
701
702 def ProcessReceptorOptions():
703 """Process receptor options. """
704
705 ReceptorFile, ReceptorMapsPrefix, UseReceptorFile, UseReceptorMapsPrefix = [None, None, False, False]
706
707 Receptor = Options["--receptor"]
708 if os.path.isfile(Receptor):
709 ReceptorFile = Receptor
710 UseReceptorFile = True
711 else:
712 ReceptorMapsPrefix = Receptor
713 UseReceptorMapsPrefix = True
714
715 OptionsInfo["Receptor"] = Receptor
716
717 OptionsInfo["ReceptorFile"] = ReceptorFile
718 OptionsInfo["UseReceptorFile"] = UseReceptorFile
719
720 OptionsInfo["ReceptorMapsPrefix"] = ReceptorMapsPrefix
721 OptionsInfo["UseReceptorMapsPrefix"] = UseReceptorMapsPrefix
722
723 UseReceptorFlexFile = False
724 ReceptorFlexFile = None
725 if not re.match("^None$", Options["--receptorFlexFile"],re.I):
726 UseReceptorFlexFile = True
727 ReceptorFlexFile = Options["--receptorFlexFile"]
728 OptionsInfo["ReceptorFlexFile"] = ReceptorFlexFile
729 OptionsInfo["UseReceptorFlexFile"] = UseReceptorFlexFile
730
731 if OptionsInfo["UseAD4Forcefield"]:
732 if OptionsInfo["UseReceptorFile"]:
733 MiscUtil.PrintError("The value specified, %s, for option \"-r, --receptor\" is not valid for %s forcefield. Supported value: Affinity maps prefix." % (OptionsInfo["Receptor"], Options["--forcefield"]))
734
735 def ProcessForcefieldWeightParamatersOption(ParamsOptionName, ParamsOptionValue, ParamsDefaultInfo = None):
736 """Process forcefield weight paramaters option."""
737
738 ParamsInfo = {"AD4Vdw": 0.1662, "AD4HydrogenBond": 0.1209, "AD4Electrostatic": 0.1406, "AD4Desolvation": 0.1322, "AD4GlueLinearAttraction": 50.0, "AD4Rot": 0.2983,
739 "VinaGaussian1": -0.035579, "VinaGaussian2": -0.005156, "VinaRepulsion": 0.840245, "VinaHydrophobic": -0.035069, "VinaHydrogenBond": -0.587439, "VinaGlueLinearAttraction": 50.0, "VinaRot": 0.05846,
740 "VinardoGaussian1": -0.045, "VinardoRepulsion": 0.8, "VinardoHydrophobic": -0.035, "VinardoHydrogenBond": -0.600, "VinardoGlueLinearAttraction": 50.0, "VinardoRot": 0.05846}
741
742 # Setup a canonical paramater names...
743 ValidParamNames = []
744 CanonicalParamNamesMap = {}
745 for ParamName in sorted(ParamsInfo):
746 ValidParamNames.append(ParamName)
747 CanonicalParamNamesMap[ParamName.lower()] = ParamName
748
749 # Update default values...
750 if ParamsDefaultInfo is not None:
751 for ParamName in ParamsDefaultInfo:
752 if ParamName not in ParamsInfo:
753 MiscUtil.PrintError("The default parameter name, %s, specified using \"%s\" to function ProcessForcefieldWeightParamatersOption is not a valid name. Supported parameter names: %s" % (ParamName, ParamsDefaultInfo, " ".join(ValidParamNames)))
754 ParamsInfo[ParamName] = ParamsDefaultInfo[ParamName]
755
756 if re.match("^auto$", ParamsOptionValue, re.I):
757 return ParamsInfo
758
759 ParamsOptionValueWords = ParamsOptionValue.split(",")
760 if len(ParamsOptionValueWords) % 2:
761 Miscutil.PrintError("The number of comma delimited paramater names and values, %d, specified using \"%s\" option must be an even number." % (len(ParamsOptionValueWords), ParamsOptionName))
762
763 # Validate paramater name and value pairs...
764 for Index in range(0, len(ParamsOptionValueWords), 2):
765 Name = ParamsOptionValueWords[Index].strip()
766 Value = ParamsOptionValueWords[Index + 1].strip()
767
768 CanonicalName = Name.lower()
769 if CanonicalName not in CanonicalParamNamesMap:
770 MiscUtil.PrintError("The parameter name, %s, specified using \"%s\" is not a valid name. Supported parameter names: %s" % (Name, ParamsOptionName, " ".join(ValidParamNames)))
771
772 ParamName = CanonicalParamNamesMap[CanonicalName]
773
774 if not MiscUtil.IsFloat(Value):
775 MiscUtil.PrintError("The parameter value, %s, specified for parameter name, %s, using \"%s\" must be a float." % (Value, Name, ParamsOptionName))
776 ParamValue = float(Value)
777
778 ParamsInfo[ParamName] = ParamValue
779
780 return ParamsInfo
781
782 def ProcessModeOption():
783 """Process mode option."""
784
785 Mode = Options["--mode"]
786 DockMode, LocalOptimizationOnlyMode, ScoreOnlyMode = [False] * 3
787 if re.match("^Dock$", Mode, re.I):
788 Mode = "Dock"
789 DockMode = True
790 elif re.match("^LocalOptimizationOnly$", Mode, re.I):
791 Mode = "LocalOptimizationOnly"
792 LocalOptimizationOnlyMode = True
793 elif re.match("^ScoreOnly$", Mode, re.I):
794 Mode = "ScoreOnly"
795 ScoreOnlyMode = True
796 else:
797 MiscUtil.PrintError("The value specified, %s, for option \"-m, --mode\" is not valid. Supported values: Dock LocalOptimizationOnly ScoreOnly" % OptionsInfo["Mode"])
798
799 OptionsInfo["Mode"] = Mode
800 OptionsInfo["DockMode"] = DockMode
801 OptionsInfo["LocalOptimizationOnlyMode"] = LocalOptimizationOnlyMode
802 OptionsInfo["ScoreOnlyMode"] = ScoreOnlyMode
803
804 def ProcessEnergyLabelOptions():
805 """Process energy label options."""
806
807 Forcefield = OptionsInfo["Forcefield"]
808
809 EnergyLabel = Options["--energyLabel"]
810 if re.match("^auto$", EnergyLabel, re.I):
811 EnergyLabel = "%s_Total_Energy (kcal/mol)" % Forcefield
812 OptionsInfo["EnergyLabel"] = EnergyLabel
813
814 EnergyLabelsList = []
815 DefaultLabels = ["%s_Intermolecular_Energy (kcal/mol)" % (Forcefield), "%s_Internal_Energy (kcal/mol)" % (Forcefield), "%s_Torsions_Energy (kcal/mol)" % (Forcefield)]
816
817 EnergyLabels = Options["--energyComponentsLabels"]
818 if not re.match("^auto$", EnergyLabels, re.I):
819 EnergyLabelsWords = EnergyLabels.split(",")
820 if len(EnergyLabelsWords) != 3:
821 MiscUtil.PrintError("The specified value, %s, for option \"--energyComponentsLabels \" is not valid. It must contain 3 text values separated by comma." % EnergyLabels)
822
823 for LabelIndex, Label in enumerate(EnergyLabelsWords):
824 Label = Label.strip()
825 if re.match("^auto$", Label, re.I):
826 Label = DefaultLabels[LabelIndex]
827
828 EnergyLabelsList.append(Label)
829 else:
830 EnergyLabelsList = DefaultLabels
831
832 OptionsInfo["EnergyComponentsLabels"] = EnergyLabels
833 OptionsInfo["EnergyComponentsLabelsList"] = EnergyLabelsList
834
835 SetupEnergyLabelsAndIndices()
836
837 def SetupEnergyLabelsAndIndices():
838 """Setup energy data labels and indices for writing out energy values. """
839
840 EnergyLabelsList = []
841 EnergyValueIndicesList = []
842
843 # Total energy is always at index 0 in the list retured by vina.score (),
844 # vina.optimize(), and vina.energies() during ScoreOnly, LocalOptimizationOnly
845 # and Dock modes...
846 EnergyLabelsList.append(OptionsInfo["EnergyLabel"])
847 EnergyValueIndicesList.append(0)
848
849 OptionsInfo["EnergyLabelsList"] = EnergyLabelsList
850 OptionsInfo["EnergyValueIndicesList"] = EnergyValueIndicesList
851
852 if not OptionsInfo["EnergyComponents"]:
853 return
854
855 # Setup energy labels and indices for energy components...
856 if OptionsInfo["ScoreOnlyMode"]:
857 # vina.score returns a list containing the following values:
858 #
859 # Vina/Vinardo FF: columns=[total, lig_inter, flex_inter, other_inter, flex_intra, lig_intra, torsions, lig_intra best pose]
860 # AutoDock FF: [total, lig_inter, flex_inter, other_inter, flex_intra, lig_intra, torsions, -lig_intra]
861 #
862 IntermolecularEnergyIndex, IntramolecularEnergyIndex, TorsionEnergyIndex = [1, 5, 6]
863 elif OptionsInfo["LocalOptimizationOnlyMode"]:
864 # vina.optimize returns a list containing the following values:
865 #
866 # Vina/Vinardo FF: [total, lig_inter, flex_inter, other_inter, flex_intra, lig_intra, torsions, lig_intra best pose]
867 # AutoDock FF: [total, lig_inter, flex_inter, other_inter, flex_intra, lig_intra, torsions, -lig_intra]
868 #
869 IntermolecularEnergyIndex, IntramolecularEnergyIndex, TorsionEnergyIndex = [1, 5, 6]
870 elif OptionsInfo["DockMode"]:
871 # vina.energies returns a list containing the following values:
872 #
873 # Vina/Vinardo FF: [total, inter, intra, torsions, intra best pose]
874 # AutoDock FF: [total, inter, intra, torsions, -intra]
875 #
876 IntermolecularEnergyIndex, IntramolecularEnergyIndex, TorsionEnergyIndex = [1, 2, 3]
877 else:
878 MiscUtil.PrintError("The value specified, %s, for option \"-m, --mode\" is not valid. Supported values: Dock LocalOptimizationOnly ScoreOnly" % OptionsInfo["Mode"])
879
880 OptionsInfo["EnergyLabelsList"].extend(OptionsInfo["EnergyComponentsLabelsList"])
881 OptionsInfo["EnergyValueIndicesList"].extend([IntermolecularEnergyIndex, IntramolecularEnergyIndex, TorsionEnergyIndex])
882
883 def ProcessOptions():
884 """Process and validate command line arguments and options."""
885
886 MiscUtil.PrintInfo("Processing options...")
887
888 # Validate options...
889 ValidateOptions()
890
891 OptionsInfo["Infile"] = Options["--infile"]
892 ParamsDefaultInfoOverride = {"RemoveHydrogens": False}
893 OptionsInfo["InfileParams"] = MiscUtil.ProcessOptionInfileParameters("--infileParams", Options["--infileParams"], InfileName = Options["--infile"], ParamsDefaultInfo = ParamsDefaultInfoOverride)
894
895 OptionsInfo["Outfile"] = Options["--outfile"]
896 OptionsInfo["OutfileParams"] = MiscUtil.ProcessOptionOutfileParameters("--outfileParams", Options["--outfileParams"])
897
898 FileDir, FileName, FileExt = MiscUtil.ParseFileName(Options["--outfile"])
899 OptionsInfo["OutfileFlexRes"] = "%s_Flex_Receptor.%s" % (FileName, FileExt)
900
901 GridCenterLigandFile, GridCenterList, GridCenterByLigandFile= [None, None, False]
902 GridCenter = Options["--gridCenter"]
903 if re.search(",", GridCenter, re.I):
904 GridCenterList = [float(Value) for Value in GridCenter.split(",")]
905 else:
906 GridCenterByLigandFile = True
907 GridCenterLigandFile = GridCenter
908 GridCenterList = ComputeGridCenter(GridCenterLigandFile)
909
910 OptionsInfo["GridCenter"] = GridCenter
911 OptionsInfo["GridCenterList"] = GridCenterList
912 OptionsInfo["GridCenterLigandFile"] = GridCenterLigandFile
913 OptionsInfo["GridCenterByLigandFile"] = GridCenterByLigandFile
914
915 OptionsInfo["EnergyComponents"] = True if re.match("^yes$", Options["--energyComponents"], re.I) else False
916
917 OptionsInfo["EnergyRange"] = float(Options["--energyRange"])
918 OptionsInfo["Exhaustiveness"] = int(Options["--exhaustiveness"])
919
920 Forcefield = Options["--forcefield"]
921 UseAD4Forcefield, UseVinaForcefield, UseVinardoForcefield = [False, False, False]
922 if re.match("^AD4$", Forcefield, re.I):
923 Forcefield = "AD4"
924 UseAD4Forcefield = True
925 elif re.match("^Vina$", Forcefield, re.I):
926 Forcefield = "Vina"
927 UseVinaForcefield = True
928 elif re.match("^Vinardo$", Forcefield, re.I):
929 Forcefield = "Vinardo"
930 UseVinardoForcefield = True
931 else:
932 MiscUtil.PrintError("The value specified, %s, for option \"-f, --forcefield\" is not valid. Supported values: AD4 Vina Vinardo")
933 OptionsInfo["Forcefield"] = Forcefield
934 OptionsInfo["UseAD4Forcefield"] = UseAD4Forcefield
935 OptionsInfo["UseVinaForcefield"] = UseVinaForcefield
936 OptionsInfo["UseVinardoForcefield"] = UseVinardoForcefield
937
938 OptionsInfo["ForcefieldWeightParams"] = ProcessForcefieldWeightParamatersOption('--forcefieldWeightParams', Options["--forcefieldWeightParams"])
939 OptionsInfo["ForcefieldWeightParamsSpecified"] = False if re.match("^auto$", Options["--forcefieldWeightParams"], re.I) else True
940
941 GridSize = Options["--gridSize"]
942 GridSizeList = [float(Value) for Value in GridSize.split(",")]
943 OptionsInfo["GridSize"] = GridSize
944 OptionsInfo["GridSizeList"] = GridSizeList
945
946 OptionsInfo["GridSpacing"] = float(Options["--gridSpacing"])
947
948 OptionsInfo["MaxEvaluations"] = int(Options["--maxEvaluations"])
949 OptionsInfo["MinRMSD"] = float(Options["--minRMSD"])
950
951 MergeHydrogens = Options["--mergeHydrogens"]
952 if re.match("^auto$", MergeHydrogens, re.I):
953 MergeHydrogens = True if OptionsInfo["UseAD4Forcefield"] else False
954 else:
955 MergeHydrogens = True if re.match("^yes$", MergeHydrogens, re.I) else False
956 OptionsInfo["MergeHydrogens"] = MergeHydrogens
957
958 ProcessModeOption()
959
960 OptionsInfo["MPMode"] = True if re.match("^yes$", Options["--mp"], re.I) else False
961 OptionsInfo["MPParams"] = MiscUtil.ProcessOptionMultiprocessingParameters("--mpParams", Options["--mpParams"])
962
963 OptionsInfo["NumPoses"] = int(Options["--numPoses"])
964
965 if re.match("^auto$", Options["--numThreads"], re.I):
966 NumThreads = 1 if OptionsInfo["MPMode"] else 0
967 else:
968 NumThreads = int(Options["--numThreads"])
969 OptionsInfo["NumThreads"] = NumThreads
970
971 OptionsInfo["Precision"] = int(Options["--precision"])
972 OptionsInfo["RandomSeed"] = int(Options["--randomSeed"])
973
974 OptionsInfo["SkipRefinement"] = True if re.match("^yes$", Options["--skipRefinement"], re.I) else False
975
976 OptionsInfo["ValidateMolecules"] = True if re.match("^yes$", Options["--validateMolecules"], re.I) else False
977
978 if re.match("^auto$", Options["--vinaVerbosity"], re.I):
979 VinaVerbosity = 0 if OptionsInfo["MPMode"] else 1
980 else:
981 VinaVerbosity = int(Options["--vinaVerbosity"])
982 OptionsInfo["VinaVerbosity"] = VinaVerbosity
983
984 OptionsInfo["Overwrite"] = Options["--overwrite"]
985 OptionsInfo["QuietMode"] = True if re.match("^yes$", Options["--quiet"], re.I) else False
986
987 ProcessEnergyLabelOptions()
988 ProcessReceptorOptions()
989
990 def RetrieveOptions():
991 """Retrieve command line arguments and options."""
992
993 # Get options...
994 global Options
995 Options = docopt(_docoptUsage_)
996
997 # Set current working directory to the specified directory...
998 WorkingDir = Options["--workingdir"]
999 if WorkingDir:
1000 os.chdir(WorkingDir)
1001
1002 # Handle examples option...
1003 if "--examples" in Options and Options["--examples"]:
1004 MiscUtil.PrintInfo(MiscUtil.GetExamplesTextFromDocOptText(_docoptUsage_))
1005 sys.exit(0)
1006
1007 def ValidateOptions():
1008 """Validate option values."""
1009
1010 MiscUtil.ValidateOptionFilePath("-i, --infile", Options["--infile"])
1011 MiscUtil.ValidateOptionFileExt("-i, --infile", Options["--infile"], "sdf sd mol")
1012
1013 MiscUtil.ValidateOptionFileExt("-o, --outfile", Options["--outfile"], "sdf sd")
1014 MiscUtil.ValidateOptionsOutputFileOverwrite("-o, --outfile", Options["--outfile"], "--overwrite", Options["--overwrite"])
1015 MiscUtil.ValidateOptionsDistinctFileNames("-i, --infile", Options["--infile"], "-o, --outfile", Options["--outfile"])
1016
1017 FileDir, FileName, FileExt = MiscUtil.ParseFileName(Options["--receptor"])
1018 if not MiscUtil.IsEmpty(FileExt):
1019 MiscUtil.ValidateOptionFilePath("-r, --receptor", Options["--receptor"])
1020 MiscUtil.ValidateOptionFileExt("-r, --receptor", Options["--receptor"], "pdbqt")
1021 else:
1022 AffinityMapFiles = glob.glob("%s.*.map" % Options["--receptor"])
1023 if len(AffinityMapFiles) == 0:
1024 MiscUtil.PrintError("The receptor affinity map files, %s.*.map, corresponding to maps prefix, %s, specified using option, \"-r, --receptor\" option don't exist." % (Options["--receptor"], Options["--receptor"]))
1025
1026 if not re.match("^None$", Options["--receptorFlexFile"],re.I):
1027 MiscUtil.ValidateOptionFilePath("--receptorFlexFile", Options["--receptorFlexFile"])
1028 MiscUtil.ValidateOptionFileExt("--receptorFlexFile", Options["--receptorFlexFile"], "pdbqt")
1029 if os.path.isfile(Options["--receptor"]):
1030 MiscUtil.ValidateOptionsDistinctFileNames("-r, --receptor", Options["--receptor"], "--receptorFlexFile", Options["--receptorFlexFile"])
1031
1032 if re.search(",", Options["--gridCenter"], re.I):
1033 MiscUtil.ValidateOptionNumberValues("-g, --gridCenter", Options["--gridCenter"], 3, ",", "float", {">=" : 0.0})
1034 else:
1035 MiscUtil.ValidateOptionFilePath("-g, --gridCenter", Options["--gridCenter"])
1036 MiscUtil.ValidateOptionFileExt("-g, --gridCenter", Options["--gridCenter"], "sdf sd mol pdb")
1037
1038 MiscUtil.ValidateOptionTextValue("--energyComponents", Options["--energyComponents"], "yes no")
1039
1040 MiscUtil.ValidateOptionFloatValue("--energyRange", Options["--energyRange"], {">": 0})
1041
1042 MiscUtil.ValidateOptionIntegerValue("--exhaustiveness", Options["--exhaustiveness"], {">": 0})
1043
1044 MiscUtil.ValidateOptionTextValue("-f, --forcefield", Options["--forcefield"], "AD4 Vina Vinardo")
1045
1046 MiscUtil.ValidateOptionNumberValues("--gridSize", Options["--gridSize"], 3, ",", "float", {">" : 0.0})
1047 MiscUtil.ValidateOptionFloatValue("--gridSpacing", Options["--gridSpacing"], {">": 0})
1048
1049 MiscUtil.ValidateOptionIntegerValue("--maxEvaluations", Options["--maxEvaluations"], {">=": 0})
1050 MiscUtil.ValidateOptionFloatValue("--minRMSD", Options["--minRMSD"], {">": 0})
1051
1052 MiscUtil.ValidateOptionTextValue("--mergeHydrogens", Options["--mergeHydrogens"], "yes no auto")
1053
1054 MiscUtil.ValidateOptionTextValue("-m, --mode", Options["--mode"], "Dock LocalOptimizationOnly ScoreOnly")
1055
1056 MiscUtil.ValidateOptionTextValue("--mp", Options["--mp"], "yes no")
1057
1058 MiscUtil.ValidateOptionIntegerValue("--numPoses", Options["--numPoses"], {">": 0})
1059
1060 if not re.match("^auto$", Options["--numThreads"], re.I):
1061 MiscUtil.ValidateOptionIntegerValue("--numThreads", Options["--numThreads"], {">": 0})
1062
1063 MiscUtil.ValidateOptionIntegerValue("-p, --precision", Options["--precision"], {">": 0})
1064 MiscUtil.ValidateOptionIntegerValue("--randomSeed", Options["--randomSeed"], {})
1065
1066 MiscUtil.ValidateOptionTextValue("--skipRefinement", Options["--skipRefinement"], "yes no")
1067
1068 MiscUtil.ValidateOptionTextValue("-v, --validateMolecules", Options["--validateMolecules"], "yes no")
1069
1070 if not re.match("^auto$", Options["--vinaVerbosity"], re.I):
1071 MiscUtil.ValidateOptionIntegerValue("--vinaVerbosity", Options["--vinaVerbosity"], {">=": 0})
1072 MiscUtil.ValidateOptionTextValue("--vinaVerbosity", Options["--vinaVerbosity"], "0 1 2")
1073
1074 MiscUtil.ValidateOptionTextValue("-q, --quiet", Options["--quiet"], "yes no")
1075
1076 # Setup a usage string for docopt...
1077 _docoptUsage_ = """
1078 VinaPerformDocking.py - Perform docking
1079
1080 Usage:
1081 VinaPerformDocking.py [--energyComponents <yes or no>] [--energyComponentsLabels <Label1, label2, Label3>]
1082 [--energyLabel <text>] [--energyRange <number>] [--exhaustiveness <number>]
1083 [--forcefield <AD4, Vina, or Vinardo>] [--forcefieldWeightParams <Name,Value,...>] [--gridSpacing <number>]
1084 [--gridSize <xsize,ysize,zsize>] [--infileParams <Name,Value,...>] [--maxEvaluations <number>]
1085 [--mergeHydrogens <yes or no>] [--minRMSD <number>] [--mode <Dock, LocalOptimizationOnly, ScoreOnly>] [--mp <yes or no>]
1086 [--mpParams <Name, Value,...>] [--numPoses <number>] [--numThreads <number>] [--outfileParams <Name,Value,...> ]
1087 [--overwrite] [--precision <number>] [--quiet <yes or no>] [--randomSeed <number>] [--receptorFlexFile <receptor flex file>]
1088 [--skipRefinement <yes or no>] [--validateMolecules <yes or no>] [--vinaVerbosity <number>]
1089 [-w <dir>] -g <RefLigandFile or x,y,z> -r <receptorfile or maps prerfix> -i <infile> -o <outfile>
1090 VinaPerformDocking.py -h | --help | -e | --examples
1091
1092 Description:
1093 Dock molecules against a protein receptor using AutoDock Vina [Ref 168, 169].
1094 The molecules must have 3D coordinates in input file. In addition, the hydrogens
1095 must be present for all molecules in input file.
1096
1097 No protein receptor preparation is performed during docking. It must be prepared
1098 employing standalone scripts available as part of AutoDock Vina. You may
1099 optionally specify flexible residues in the binding pocket to prepare a flexible
1100 receptor file and employ it for docking molecules along with the fixed receptor
1101 file.
1102
1103 The following three forecefileds are available to score molecules: AD4
1104 (AutoDock4), Vina, and Vinardo (Vina RaDii Optimized) [Ref 170].
1105
1106 The supported input file formats are shown below:
1107
1108 Rigid/Flexible protein receptor files - PDBQT(.pdbqt)
1109 Reference ligand file: PDB(.pdb), Mol (.mol), SD (.sdf, .sd)
1110
1111 Input molecules file - Mol (.mol), SD (.sdf, .sd)
1112
1113 The supported output file format is: SD (.sdf, .sd).
1114
1115 The following output files are generated:
1116
1117 <OutfileRoot>.<OutfileExt> - Docked/scored molecules
1118 <OutfileRoot>_Flex_Receptor.<OutfileExt> - Docked poses for flexible
1119 residues
1120
1121 The flexible receptor output file contains docked poses corresponding to
1122 flexible residues. It is only generated during 'Dock' value of '-m, --mode'
1123 option. The number of poses in this file matches those written to the
1124 output file containing docked molecules.
1125
1126 Options:
1127 --energyComponents <yes or no> [default: no]
1128 Write out binding energy components of the total binding energy docking
1129 score to outfile. The following three energy components are written to
1130 outfile: intermolecula energy, internal energy, and torsions energy.
1131 --energyComponentsLabels <Label1, label2, Label3> [default: auto]
1132 A triplet of comma delimited values corresponding to energy data field
1133 labels for writing out the binding energy components to outfile. You must
1134 specify all three values. A value of 'None' implies the use of the default
1135 labels as shown below:
1136
1137 Label1: <ForcefieldName>_Intermolecular_Energy (kcal/mol)
1138 Label2: <ForcefieldName>_Internal_Energy (kcal/mol)
1139 Label3: <ForcefieldName>_Torsions_Energy (kcal/mol)
1140
1141 --energyLabel <text> [default: auto]
1142 Energy data field label for writing out binding energy docking score to
1143 output file. Default: <ForcefieldName>_Total_Energy (kcal/mol).
1144 --energyRange <number> [default: 3.0]
1145 Maximum energy difference from the best pose during the generation of
1146 poses. Units: kcal/mol.
1147 -e, --examples
1148 Print examples.
1149 --exhaustiveness <number> [default: 8]
1150 Exhaustiveness of global MC search. The higher values make the search
1151 more exhaustive and it takes longer to complete. You may want to use
1152 '16' or '32' as the value of '--exhaustiveness ' to increase the accuracy of
1153 your pose prediction.
1154 -f, --forcefield <AD4, Vina, or Vinardo> [default: Vina]
1155 Forcefield to use for scoring. Possible values: AD4 (AutoDock 4), Vina
1156 [Ref 169, 169], or Vinardo (Vina RaDii Optimized) [Ref 170].
1157
1158 You must specify affinity maps using '-r, --receptor' option during the use
1159 of 'AD4' forcefield.
1160 --forcefieldWeightParams <Name,Value,...> [default: auto]
1161 A comma delimited list of parameter name and value pairs for forcefield
1162 scoring.
1163
1164 The supported parameter names along with their default values are
1165 are shown below for different forcefields:
1166
1167 AD4 (6 weights):
1168
1169 ad4Vdw, 0.1662, ad4HydrogenBond, 0.1209, ad4Electrostatic, 0.1406,
1170 ad4Desolvation, 0.1322, ad4GlueLinearAttraction, 50.0,
1171 ad4Rot, 0.2983
1172
1173 Vina (7 weights):
1174
1175 vinaGaussian1, -0.035579, vinaGaussian2, -0.005156,
1176 vinaRepulsion, 0.840245, vinaHydrophobic, -0.035069,
1177 vinaHydrogenBond, -0.587439, vinaGlueLinearAttraction, 50.0,
1178 vinaRot, 0.05846
1179
1180 Vinardo (6 weights):
1181
1182 vinardoGaussian1, -0.045, vinardoRepulsion, 0.8,
1183 vinardoHydrophobic, -0.035, vinardoHydrogenBond, -0.600
1184 vinardoGlueLinearAttraction, 50.0, vinardoRot, 0.05846
1185
1186 The glue weight parameter corresponds to linear attraction for macrocycle
1187 closure and has the same value for AD4, Vina, and Vinardo. The rot weight
1188 has the same value for Vina and Vinardo.
1189 -g, --gridCenter <RefLigandFile or x,y,z>
1190 Reference ligand file for calculating the docking grid center or a triplet
1191 of comma delimited values in Angstrom corresponding to grid center.
1192
1193 This is required option. However, it is ignored during the specification
1194 of maps prefix for '-r, --receptor' option.
1195 --gridSize <xsize,ysize,zsize> [default: 25.0, 25.0, 25.0]
1196 Docking grid size in Angstrom.
1197 --gridSpacing <number> [default: 0.375]
1198 Docking grid spacing in Angstrom.
1199 -h, --help
1200 Print this help message.
1201 -i, --infile <infile>
1202 Input file name containing molecules for docking against a receptor. The
1203 molecules must have 3D coordinates in input file. In addition, the hydrogens
1204 must be present for all molecules in input file. The input file may contain 3D
1205 conformers.
1206 --infileParams <Name,Value,...> [default: auto]
1207 A comma delimited list of parameter name and value pairs for reading
1208 molecules from files. The supported parameter names for different file
1209 formats, along with their default values, are shown below:
1210
1211 SD, MOL: removeHydrogens,no,sanitize,yes,strictParsing,yes
1212
1213 --maxEvaluations <number> [default: 0]
1214 Maximum number of evaluations to perform for each MC run during docking.
1215 By default, its value is set 0 and the number of MC evaluations is determined
1216 using heuristic rules.
1217 --mergeHydrogens <yes or no> [default: auto]
1218 Merge hydrogens during preparation of molecules for docking. The hydrogens
1219 are automatically merged during 'AD4' value of '-f, --forcefield' option and its
1220 value is set to 'yes'. Otherwise, it's set to 'no'.
1221 --minRMSD <number> [default: 1.0]
1222 Minimum RMSD between output poses in Angstrom.
1223 -m, --mode <Dock, LocalOptimizationOnly, ScoreOnly> [default: Dock]
1224 Dock molecules or simply score molecules without performing any docking.
1225 The supported values along with a brief explanation of the expected
1226 behavior are shown below:
1227
1228 Dock: Global search along with local optimization and scoring after
1229 docking
1230 LocalOptimizationOnly: Local optimization and scoring without any docking
1231 ScoreOnly: Scoring without any local optimizatiom and docking
1232
1233 The 'ScoreOnly" allows you to score 3D moleculed from input file which
1234 are already positioned in a binding pocket of a receptor.
1235 --mp <yes or no> [default: no]
1236 Use multiprocessing.
1237
1238 By default, input data is retrieved in a lazy manner via mp.Pool.imap()
1239 function employing lazy RDKit data iterable. This allows processing of
1240 arbitrary large data sets without any additional requirements memory.
1241
1242 All input data may be optionally loaded into memory by mp.Pool.map()
1243 before starting worker processes in a process pool by setting the value
1244 of 'inputDataMode' to 'InMemory' in '--mpParams' option.
1245
1246 A word to the wise: The default 'chunkSize' value of 1 during 'Lazy' input
1247 data mode may adversely impact the performance. The '--mpParams' section
1248 provides additional information to tune the value of 'chunkSize'.
1249 --mpParams <Name,Value,...> [default: auto]
1250 A comma delimited list of parameter name and value pairs to configure
1251 multiprocessing.
1252
1253 The supported parameter names along with their default and possible
1254 values are shown below:
1255
1256 chunkSize, auto
1257 inputDataMode, Lazy [ Possible values: InMemory or Lazy ]
1258 numProcesses, auto [ Default: mp.cpu_count() ]
1259
1260 These parameters are used by the following functions to configure and
1261 control the behavior of multiprocessing: mp.Pool(), mp.Pool.map(), and
1262 mp.Pool.imap().
1263
1264 The chunkSize determines chunks of input data passed to each worker
1265 process in a process pool by mp.Pool.map() and mp.Pool.imap() functions.
1266 The default value of chunkSize is dependent on the value of 'inputDataMode'.
1267
1268 The mp.Pool.map() function, invoked during 'InMemory' input data mode,
1269 automatically converts RDKit data iterable into a list, loads all data into
1270 memory, and calculates the default chunkSize using the following method
1271 as shown in its code:
1272
1273 chunkSize, extra = divmod(len(dataIterable), len(numProcesses) * 4)
1274 if extra: chunkSize += 1
1275
1276 For example, the default chunkSize will be 7 for a pool of 4 worker processes
1277 and 100 data items.
1278
1279 The mp.Pool.imap() function, invoked during 'Lazy' input data mode, employs
1280 'lazy' RDKit data iterable to retrieve data as needed, without loading all the
1281 data into memory. Consequently, the size of input data is not known a priori.
1282 It's not possible to estimate an optimal value for the chunkSize. The default
1283 chunkSize is set to 1.
1284
1285 The default value for the chunkSize during 'Lazy' data mode may adversely
1286 impact the performance due to the overhead associated with exchanging
1287 small chunks of data. It is generally a good idea to explicitly set chunkSize to
1288 a larger value during 'Lazy' input data mode, based on the size of your input
1289 data and number of processes in the process pool.
1290
1291 The mp.Pool.map() function waits for all worker processes to process all
1292 the data and return the results. The mp.Pool.imap() function, however,
1293 returns the the results obtained from worker processes as soon as the
1294 results become available for specified chunks of data.
1295
1296 The order of data in the results returned by both mp.Pool.map() and
1297 mp.Pool.imap() functions always corresponds to the input data.
1298 -n, --numPoses <number> [default: 1]
1299 Number of docked poses to generate for each molecule and write out
1300 to output file. This option is only valid for 'Dock' value of '-m, --mode'
1301 option.
1302 --numThreads <number> [default: auto]
1303 Number of threads/CPUs to use for MC search calculation in Vina. The
1304 default value is set to 1 during multiprocessing for 'Yes' value of '--mp'
1305 option. Otherwise, it's set to 0 and rely on Vina detect and use all
1306 available CPUs for multi-threading.
1307 -o, --outfile <outfile>
1308 Output file name for writing out molecules. The flexible receptor residues
1309 are written to <OutfileRoot>_Flex_Receptor.<OutfileExt>.
1310 --outfileParams <Name,Value,...> [default: auto]
1311 A comma delimited list of parameter name and value pairs for writing
1312 molecules to files. The supported parameter names for different file
1313 formats, along with their default values, are shown below:
1314
1315 SD: kekulize,yes,forceV3000,no
1316
1317 --overwrite
1318 Overwrite existing files.
1319 --precision <number> [default: 2]
1320 Floating point precision for writing energy values.
1321 -q, --quiet <yes or no> [default: no]
1322 Use quiet mode. The warning and information messages will not be printed.
1323 --randomSeed <number> [default: 0]
1324 Random seed for MC search calculations. A value of zero implies it's randomly
1325 chosen by Vina during the calculation.
1326 -r, --receptor <receptor file or maps prerfix>
1327 Protein receptor file name or prefix for affinity map files corresponding
1328 to the fixed portion of the receptor.
1329
1330 You must specify affinity map files for 'AD4' forcefield. The affinity map
1331 files correspond to <MapsPrefix>.*.map and must be present
1332
1333 The supported receptor file format is PDBQT (.pdbqt). It must contain a
1334 prepared protein receptor ready for docking. You may prepare a PDBQT
1335 receptor file from a PDB file employing the command line scripts available
1336 with AutoDock Vina and Meeko. For example: prepare_receptor,
1337 prepare_flexreceptor.py, or or mk_make_recepror.py.
1338
1339 You may want to perform the following steps to clean up your PDB file
1340 before generating a PDBQT receptor file: Remove extraneous molecules
1341 such as solvents, ions, and ligand etc.; Extract data for a chain containing
1342 the binding pocket of interest.
1343 --receptorFlexFile <receptor flex file> [default: none]
1344 Protein receptor file name corresponding to the flexible portion of the
1345 receptor. The supported receptor file format is PDBQT (.pdbqt). It must
1346 contain a prepared protein receptor ready for docking. You may prepare
1347 a flexible PDBQT receptor file from a PDB file employing the command line
1348 script prepare_flexreceptor or mk_make_recepror.py available with Autodock
1349 Vina and Meeko.
1350 --skipRefinement <yes or no> [default: no]
1351 Skip refinement. Vina is initialized to skip the use of explicit receptor atoms,
1352 instead of precalculated grids, during the following three docking modes:
1353 Dock, LocalOptimizationOnly, and ScoreOnly.
1354 -v, --validateMolecules <yes or no> [default: yes]
1355 Validate molecules for docking. The input molecules must have 3D coordinates
1356 and the hydrogens must be present. You may skip validation of molecules in
1357 input file containing all valid molecules.
1358 --vinaVerbosity <number> [default: auto]
1359 Verbosity level for running Vina. Possible values: 0 - No output; 1 - Normal
1360 output; 2 - Verbose. Default: 0 during multiprocessing for 'Yes' value of '--mp'
1361 option; otherwise, its set to 1. A non-zero value is not recommended during
1362 multiprocessing. It doesn't work due to the mingling of the Vina output
1363 from multiple processes.
1364 -w, --workingdir <dir>
1365 Location of working directory which defaults to the current directory.
1366
1367 Examples:
1368 To dock molecules using Vina forcefield against a prepared receptor
1369 corresponding to chain A extracted from PDB file 1R4L.pdb, multi-threading
1370 across all available CPUs during Vina MC search, calculating grid center form
1371 a reference ligand file, using grid box size of 25 Angstrom, generating one
1372 pose for each molecule, and write out a SD file containing docked molecules,
1373 type:
1374
1375 % VinaPerformDocking.py -g SampleACE2RefLigand.pdb
1376 -r SampleACE2Receptor.pdbqt -i SampleACE2Ligands.sdf
1377 -o SampleACE2LigandsOut.sdf
1378
1379 To run the first example for generating multiple docked poses for each
1380 molecule and write out all energy terms to a SD file, type:
1381
1382 % VinaPerformDocking.py -g SampleACE2RefLigand.pdb
1383 -r SampleACE2Receptor.pdbqt --numPoses 5 --energyComponents yes
1384 -i SampleACE2Ligands.sdf -o SampleACE2LigandsOut.sdf
1385
1386 To run the first example for docking molecules using Vinardo forcefield and write
1387 out a SD file, type:
1388
1389 % VinaPerformDocking.py -f Vinardo -g SampleACE2RefLigand.pdb
1390 -r SampleACE2Receptor.pdbqt -i SampleACE2Ligands.sdf
1391 -o SampleACE2LigandsOut.sdf
1392
1393 To run the first example for docking molecules using AD4 forcefield relying on
1394 the presence of affinity maps in the working directory and write out a SD file,
1395 type:
1396
1397 % VinaPerformDocking.py -f AD4 -g SampleACE2RefLigand.pdb
1398 -r SampleACE2Receptor -i SampleACE2Ligands.sdf
1399 -o SampleACE2LigandsOut.sdf
1400
1401 To run the first example for docking molecules using a set of explicit values
1402 for grid dimensions and write out a SD file, type:
1403
1404 % VinaPerformDocking.py -g "41.399, 5.851, 28.082"
1405 --gridSize "25.0, 25.0, 25.0" --gridSpacing 0.375
1406 -r SampleACE2Receptor.pdbqt -i SampleACE2Ligands.sdf
1407 -o SampleACE2LigandsOut.sdf
1408
1409 To run the first example for only scoring molecules already positioned in the
1410 binding pocket and write out a SD file, type:
1411
1412 % VinaPerformDocking.py -m ScoreOnly -g SampleACE2RefLigand.sdf
1413 -r SampleACE2Receptor.pdbqt -i SampleACE2RefLigandWithHs.sdf
1414 -o SampleACE2LigandsOut.sdf
1415
1416 To run the first example for docking molecules to increase the accuracy of
1417 pose predictions and write out a SD file, type:
1418
1419 % VinaPerformDocking.py -g SampleACE2RefLigand.pdb
1420 -r SampleACE2Receptor.pdbqt --exhaustiveness 24
1421 -i SampleACE2Ligands.sdf -o SampleACE2LigandsOut.sdf
1422
1423 To run the first example in multiprocessing mode on all available CPUs
1424 without loading all data into memory, a single thread for Vina docking, and
1425 write out a SD file, type:
1426
1427 % VinaPerformDocking.py -g SampleACE2RefLigand.pdb
1428 -r SampleACE2Receptor.pdbqt --mp yes
1429 -i SampleACE2Ligands.sdf -o SampleACE2LigandsOut.sdf
1430
1431 To run the first example in multiprocessing mode on all available CPUs
1432 by loading all data into memory, a single thread for Vina, and write out
1433 a SD file, type:
1434
1435 % VinaPerformDocking.py -g SampleACE2RefLigand.pdb
1436 -r SampleACE2Receptor.pdbqt --mp yes --mpParams "inputDataMode,
1437 InMemory" -i SampleACE2Ligands.sdf -o SampleACE2LigandsOut.sdf
1438
1439 To run the first example in multiprocessing mode on specific number of CPUs
1440 and chunk size without loading all data into memory along with a specific number
1441 of threads for Vina docking and write out a SD file, type:
1442
1443 % VinaPerformDocking.py -g SampleACE2RefLigand.pdb
1444 -r SampleACE2Receptor.pdbqt --mp yes --mpParams "inputDataMode,lazy,
1445 numProcesses,4,chunkSize,2" --numThreads 2
1446 -i SampleACE2Ligands.sdf -o SampleACE2LigandsOut.sdf
1447
1448 To run the first example for docking molecules employing a flexible portion
1449 of the receptor corresponding to ARG273 and write out a SD file, type:
1450
1451 % VinaPerformDocking.py -g SampleACE2RefLigand.pdb
1452 -r SampleACE2RigidReceptor.pdbqt
1453 --receptorFlexFile SampleACE2FlexReceptor.pdbqt
1454 -i SampleACE2Ligands.sdf -o SampleACE2LigandsOut.sdf
1455
1456 To run the first example for docking molecules using specified parameters and
1457 write out a SD file, type:
1458
1459 % VinaPerformDocking.py -g "41.399, 5.851, 28.082"
1460 --gridSize "25.0, 25.0, 25.0" --gridSpacing 0.375 --energyComponents
1461 yes --exhaustiveness 32 --forcefield Vina --mode dock --numPoses 2
1462 --numThreads 4 --randomSeed 42 --validateMolecules no
1463 --vinaVerbosity 0 -r SampleACE2Receptor.pdbqt
1464 -i SampleACE2Ligands.sdf -o SampleACE2LigandsOut.sdf
1465
1466 Author:
1467 Manish Sud(msud@san.rr.com)
1468
1469 Acknowledgments:
1470 Diogo Santos-Martins and Stefano Forli
1471
1472 See also:
1473 PyMOLConvertLigandFileFormat.py, PyMOLExtractSelection.py,
1474 PyMOLInfoMacromolecules.py, PyMOLVisualizeMacromolecules.py,
1475 RDKitConvertFileFormat.py, RDKitEnumerateTautomers.py,
1476 RDKitGenerateConformers.py, RDKitPerformMinimization.py,
1477 RDKitPerformConstrainedMinimization.py
1478
1479 Copyright:
1480 Copyright (C) 2024 Manish Sud. All rights reserved.
1481
1482 The functionality available in this script is implemented using AutoDockVina
1483 and Meeko, open source software packages for docking, and RDKit, an open
1484 source toolkit for cheminformatics developed by Greg Landrum.
1485
1486 This file is part of MayaChemTools.
1487
1488 MayaChemTools is free software; you can redistribute it and/or modify it under
1489 the terms of the GNU Lesser General Public License as published by the Free
1490 Software Foundation; either version 3 of the License, or (at your option) any
1491 later version.
1492
1493 """
1494
1495 if __name__ == "__main__":
1496 main()
