Source code for ssbio.protein.structure.utils.dock

"""
DOCK
====
"""

import os
import pandas as pd
import os.path as op
import logging
from ssbio.core.object import Object
import ssbio.utils

log = logging.getLogger(__name__)


[docs]class DOCK(Object): """Class to prepare a structure file for docking with DOCK6. Attributes: """ def __init__(self, structure_id, pdb_file, amb_file, flex1_file, flex2_file, root_dir=None): """Initialize a DOCK6 project. Args: """ super(DOCK, self).__init__(id=structure_id, description='DOCK6 preparation') self._root_dir = None self.structure_path = pdb_file if root_dir: self.root_dir = root_dir else: self.root_dir = self.structure_dir self.dockprep_path = None self.receptormol2_path = None self.receptorpdb_path = None self.dms_path = None self.sphgen_path = None self.bindingsite_path = None self.sphsel_path = None self.box_path = None self.grid_path = None self.dock_flexible_outfile = None self.dock_flexible_scored_result = None self.dock_flexible_conformers_result = None self.amb_file = amb_file self.flex1_file = flex1_file self.flex2_file = flex2_file log.debug('{}: created DOCK6 project folder at {}'.format(structure_id, self.dock_dir)) @property def root_dir(self): """str: Directory where DOCK project folder is located""" return self._root_dir @root_dir.setter def root_dir(self, path): if not path: raise ValueError('No path specified') if not op.exists(path): raise ValueError('{}: folder does not exist'.format(path)) if self._root_dir: log.debug('Changing root directory of DOCK project for "{}" from {} to {}'.format(self.id, self.root_dir, path)) if not op.exists(op.join(path, self.id)): raise IOError('DOCK project "{}" does not exist in folder {}'.format(self.id, path)) self._root_dir = path for d in [self.dock_dir]: ssbio.utils.make_dir(d) @property def dock_dir(self): """str: DOCK folder""" if self.root_dir: return op.join(self.root_dir, self.id + '_DOCK') else: log.warning('Root directory not set') return None @property def structure_dir(self): if not self._structure_dir: raise OSError('No structure folder set') return self._structure_dir @structure_dir.setter def structure_dir(self, path): if path and not op.exists(path): raise OSError('{}: folder does not exist'.format(path)) self._structure_dir = path @property def structure_path(self): if not self.structure_file: raise OSError('Metadata file not loaded') path = op.join(self.structure_dir, self.structure_file) if not op.exists(path): raise OSError('{}: file does not exist'.format(path)) return path @structure_path.setter def structure_path(self, path): """Provide pointers to the paths of the structure file Args: path: Path to structure file """ if not path: self.structure_dir = None self.structure_file = None else: if not op.exists(path): raise OSError('{}: file does not exist!'.format(path)) if not op.dirname(path): self.structure_dir = '.' else: self.structure_dir = op.dirname(path) self.structure_file = op.basename(path)
[docs] def dockprep(self, force_rerun=False): """Prepare a PDB file for docking by first converting it to mol2 format. Args: force_rerun (bool): If method should be rerun even if output file exists """ log.debug('{}: running dock preparation...'.format(self.id)) prep_mol2 = op.join(self.dock_dir, '{}_prep.mol2'.format(self.id)) prep_py = op.join(self.dock_dir, "prep.py") if ssbio.utils.force_rerun(flag=force_rerun, outfile=prep_mol2): with open(prep_py, "w") as f: f.write('import chimera\n') f.write('from DockPrep import prep\n') f.write('models = chimera.openModels.list(modelTypes=[chimera.Molecule])\n') f.write('prep(models)\n') f.write('from WriteMol2 import writeMol2\n') f.write('writeMol2(models, "{}")\n'.format(prep_mol2)) cmd = 'chimera --nogui {} {}'.format(self.structure_path, prep_py) os.system(cmd) os.remove(prep_py) os.remove('{}c'.format(prep_py)) if ssbio.utils.is_non_zero_file(prep_mol2): self.dockprep_path = prep_mol2 log.debug('{}: successful dockprep execution'.format(self.dockprep_path)) else:
log.critical('{}: dockprep failed to run on PDB file'.format(self.structure_path))
[docs] def protein_only_and_noH(self, keep_ligands=None, force_rerun=False): """Isolate the receptor by stripping everything except protein and specified ligands. Args: keep_ligands (str, list): Ligand(s) to keep in PDB file force_rerun (bool): If method should be rerun even if output file exists """ log.debug('{}: running protein receptor isolation...'.format(self.id)) if not self.dockprep_path: return ValueError('Please run dockprep') receptor_mol2 = op.join(self.dock_dir, '{}_receptor.mol2'.format(self.id)) receptor_noh = op.join(self.dock_dir, '{}_receptor_noH.pdb'.format(self.id)) prly_com = op.join(self.dock_dir, "prly.com") if ssbio.utils.force_rerun(flag=force_rerun, outfile=receptor_noh): with open(prly_com, "w") as f: f.write('open {}\n'.format(self.dockprep_path)) keep_str = 'delete ~protein' if keep_ligands: keep_ligands = ssbio.utils.force_list(keep_ligands) for res in keep_ligands: keep_str += ' & ~:{} '.format(res) keep_str = keep_str.strip() + '\n' f.write(keep_str) f.write('write format mol2 0 {}\n'.format(receptor_mol2)) f.write('delete element.H\n') f.write('write format pdb 0 {}\n'.format(receptor_noh)) cmd = 'chimera --nogui {}'.format(prly_com) os.system(cmd) os.remove(prly_com) if ssbio.utils.is_non_zero_file(receptor_mol2) and ssbio.utils.is_non_zero_file(receptor_noh): self.receptormol2_path = receptor_mol2 self.receptorpdb_path = receptor_noh log.debug('{}: successful receptor isolation (mol2)'.format(self.receptormol2_path)) log.debug('{}: successful receptor isolation (pdb)'.format(self.receptorpdb_path)) else:
log.critical('{}: protein_only_and_noH failed to run on dockprep file'.format(self.dockprep_path))
[docs] def dms_maker(self, force_rerun=False): """Create surface representation (dms file) of receptor Args: force_rerun (bool): If method should be rerun even if output file exists """ log.debug('{}: running surface representation maker...'.format(self.id)) if not self.receptorpdb_path: return ValueError('Please run protein_only_and_noH') dms = op.join(self.dock_dir, '{}_receptor.dms'.format(self.id)) if ssbio.utils.force_rerun(flag=force_rerun, outfile=dms): cmd = 'dms {} -n -w 1.4 -o {}'.format(self.receptorpdb_path, dms) os.system(cmd) self.dms_path = dms if ssbio.utils.is_non_zero_file(dms): self.dms_path = dms log.debug('{}: successful dms execution'.format(self.dms_path)) else:
log.critical('{}: dms_maker failed to run on receptor file'.format(self.receptorpdb_path))
[docs] def sphgen(self, force_rerun=False): """Create sphere representation (sph file) of receptor from the surface representation Args: force_rerun (bool): If method should be rerun even if output file exists """ log.debug('{}: running sphere generation...'.format(self.id)) if not self.dms_path: return ValueError('Please run dms_maker') sph = op.join(self.dock_dir, '{}_receptor.sph'.format(self.id)) insph = op.join(self.dock_dir, 'INSPH') if ssbio.utils.force_rerun(flag=force_rerun, outfile=sph): with open(insph, "w") as f: f.write("{}\n".format(self.dms_path)) f.write("R\n") f.write("X\n") f.write("0.0\n") f.write("4.0\n") f.write("1.4\n") f.write("{}\n".format(sph)) os.chdir(self.dock_dir) cmd = "sphgen_cpp" os.system(cmd) os.remove(insph) if ssbio.utils.is_non_zero_file(sph): self.sphgen_path = sph log.debug('{}: successful sphgen execution'.format(self.sphgen_path)) else:
log.critical('{}: sphgen_cpp failed to run on dms file'.format(self.dms_path))
[docs] def binding_site_mol2(self, residues, force_rerun=False): """Create mol2 of only binding site residues from the receptor This function will take in a .pdb file (preferably the _receptor_noH.pdb file) and a string of residues (eg: '144,170,199') and delete all other residues in the .pdb file. It then saves the coordinates of the selected residues as a .mol2 file. This is necessary for Chimera to select spheres within the radius of the binding site. Args: residues (str): Comma separated string of residues (eg: '144,170,199') force_rerun (bool): If method should be rerun even if output file exists """ log.debug('{}: running binding site isolation...'.format(self.id)) if not self.receptorpdb_path: return ValueError('Please run protein_only_and_noH') prefix = self.id + '_' + 'binding_residues' mol2maker = op.join(self.dock_dir, '{}_make_mol2.py'.format(prefix)) outfile = op.join(self.dock_dir, '{}.mol2'.format(prefix)) if ssbio.utils.force_rerun(flag=force_rerun, outfile=outfile): with open(mol2maker, 'w') as mol2_maker: mol2_maker.write('#! /usr/bin/env python\n') mol2_maker.write('from chimera import runCommand\n') mol2_maker.write('runCommand("open {}")\n'.format(self.receptorpdb_path)) mol2_maker.write('runCommand("delete ~:{}")\n'.format(residues)) mol2_maker.write('runCommand("write format mol2 resnum 0 {}")\n'.format(outfile)) mol2_maker.write('runCommand("close all")') cmd = 'chimera --nogui {}'.format(mol2maker) os.system(cmd) os.remove(mol2maker) os.remove('{}c'.format(mol2maker)) if ssbio.utils.is_non_zero_file(outfile): self.bindingsite_path = outfile log.debug('{}: successful binding site isolation'.format(self.bindingsite_path)) else:
log.critical('{}: binding_site_mol2 failed to run on receptor file'.format(self.receptorpdb_path))
[docs] def sphere_selector_using_residues(self, radius, force_rerun=False): """Select spheres based on binding site residues Args: radius (int, float): Radius around binding residues to dock to force_rerun (bool): If method should be rerun even if output file exists """ log.debug('{}: running sphere selector...'.format(self.id)) if not self.sphgen_path or not self.bindingsite_path: return ValueError('Please run sphgen and binding_site_mol2') selsph = op.join(self.dock_dir, '{}_selsph_binding.sph'.format(self.id)) if ssbio.utils.force_rerun(flag=force_rerun, outfile=selsph): cmd = "sphere_selector {} {} {}".format(self.sphgen_path, self.bindingsite_path, radius) rename = "mv selected_spheres.sph {}".format(selsph) os.system(cmd) os.system(rename) if ssbio.utils.is_non_zero_file(selsph): self.sphsel_path = selsph log.debug('{}: successful sphere selection'.format(self.sphsel_path)) else:
log.critical('{}: sphere_selector_using_residues failed to run on sph file'.format(self.sphgen_path)) # def split_sph(self, force_rerun=False): # """TODO: documentation? what was this used for""" # # selsph = op.join(self.dock_dir, '{}_selsph.sph'.format(self.id)) # # if ssbio.utils.force_rerun(flag=force_rerun, outfile=selsph): # with open(self.sphgen_path, "r") as f: # text = f.read() # f.seek(0) # lines = f.readlines() # paragraphs = re.split("cluster ... number of spheres in cluster ...\n", text) # # with open(selsph, "w") as f2: # f2.write(lines[1]) # f2.write(paragraphs[1]) # # return selsph
[docs] def showbox(self, force_rerun=False): """Create the dummy PDB box around the selected spheres. Args: force_rerun (bool): If method should be rerun even if output file exists """ log.debug('{}: running box maker...'.format(self.id)) if not self.sphsel_path: return ValueError('Please run sphere_selector_using_residues') boxfile = op.join(self.dock_dir, "{}_box.pdb".format(self.id)) boxscript = op.join(self.dock_dir, "{}_box.in".format(self.id)) if ssbio.utils.force_rerun(flag=force_rerun, outfile=boxfile): with open(boxscript, "w") as f: f.write("Y\n") f.write("0\n") f.write("{}\n".format(op.basename(self.sphsel_path))) f.write("1\n") f.write("{}".format(op.basename(boxfile))) cmd = "showbox < {}".format(boxscript) os.chdir(self.dock_dir) os.system(cmd) if ssbio.utils.is_non_zero_file(boxfile): self.box_path = boxfile log.debug('{}: successful box creation'.format(self.box_path)) else:
log.critical('{}: showbox failed to run on selected spheres file'.format(self.sphsel_path))
[docs] def grid(self, force_rerun=False): """Create the scoring grid within the dummy box. Args: force_rerun (bool): If method should be rerun even if output file exists """ log.debug('{}: running grid maker...'.format(self.id)) if not self.receptormol2_path or not self.box_path: return ValueError('Please run protein_only_and_noH and showbox') gridscript = op.join(self.dock_dir, "{}_grid.in".format(self.id)) out_name = op.join(self.dock_dir, "{}_grid.out".format(self.id)) if ssbio.utils.force_rerun(flag=force_rerun, outfile=out_name): with open(gridscript, "w") as f: grid_text = """compute_grids yes grid_spacing 0.3 output_molecule no contact_score yes contact_cutoff_distance 4.5 energy_score yes energy_cutoff_distance 9999 atom_model all attractive_exponent 6 repulsive_exponent 12 distance_dielectric yes dielectric_factor 4 bump_filter yes bump_overlap 0.75 receptor_file {} box_file {} vdw_definition_file {} score_grid_prefix {}_grid """.format(op.basename(self.receptormol2_path), op.basename(self.box_path), self.amb_file, self.id) f.write(grid_text) os.chdir(self.dock_dir) cmd = "grid -i {} -o {}".format(op.basename(gridscript), op.basename(out_name)) os.system(cmd) if ssbio.utils.is_non_zero_file(out_name): self.grid_path = out_name log.debug('{}: successful grid creation'.format(self.grid_path)) else:
log.critical('{}: grid failed to run on receptor + box file'.format(self.box_path))
[docs] def do_dock6_flexible(self, ligand_path, force_rerun=False): """Dock a ligand to the protein. Args: ligand_path (str): Path to ligand (mol2 format) to dock to protein force_rerun (bool): If method should be rerun even if output file exists """ log.debug('{}: running DOCK6...'.format(self.id)) ligand_name = os.path.basename(ligand_path).split('.')[0] in_name = op.join(self.dock_dir, "{}_{}_flexdock.in".format(self.id, ligand_name)) out_name = op.join(self.dock_dir, "{}_{}_flexdock.out".format(self.id, ligand_name)) conformers_out = op.join(self.dock_dir, '{}_{}_flexdock_conformers.mol2'.format(self.id, ligand_name)) scored_out = op.join(self.dock_dir, '{}_{}_flexdock_scored.mol2'.format(self.id, ligand_name)) ranked_out = op.join(self.dock_dir, '{}_{}_flexdock_ranked.mol2'.format(self.id, ligand_name)) if ssbio.utils.force_rerun(flag=force_rerun, outfile=ranked_out): with open(in_name, "w") as f: dock_text = """ligand_atom_file {} limit_max_ligands no skip_molecule no read_mol_solvation no calculate_rmsd no use_database_filter no orient_ligand yes automated_matching yes receptor_site_file {} max_orientations 500 critical_points no chemical_matching no use_ligand_spheres no use_internal_energy yes internal_energy_rep_exp 12 flexible_ligand yes user_specified_anchor no limit_max_anchors no min_anchor_size 5 pruning_use_clustering yes pruning_max_orients 100 pruning_clustering_cutoff 100 pruning_conformer_score_cutoff 100 use_clash_overlap no write_growth_tree no bump_filter yes bump_grid_prefix {} score_molecules yes contact_score_primary no contact_score_secondary no grid_score_primary yes grid_score_secondary no grid_score_rep_rad_scale 1 grid_score_vdw_scale 1 grid_score_es_scale 1 grid_score_grid_prefix {} multigrid_score_secondary no dock3.5_score_secondary no continuous_score_secondary no descriptor_score_secondary no gbsa_zou_score_secondary no gbsa_hawkins_score_secondary no SASA_descriptor_score_secondary no amber_score_secondary no minimize_ligand yes minimize_anchor yes minimize_flexible_growth yes use_advanced_simplex_parameters no simplex_max_cycles 1 simplex_score_converge 0.1 simplex_cycle_converge 1.0 simplex_trans_step 1.0 simplex_rot_step 0.1 simplex_tors_step 10.0 simplex_anchor_max_iterations 500 simplex_grow_max_iterations 500 simplex_grow_tors_premin_iterations 0 simplex_random_seed 0 simplex_restraint_min yes simplex_coefficient_restraint 10.0 atom_model all vdw_defn_file {} flex_defn_file {} flex_drive_file {} ligand_outfile_prefix {}_{}_flexdock write_orientations no num_scored_conformers 20 write_conformations yes cluster_conformations yes rank_ligands yes """.format(ligand_path, op.basename(self.sphsel_path), op.splitext(op.basename(self.grid_path))[0], op.splitext(op.basename(self.grid_path))[0], self.amb_file, self.flex1_file, self.flex2_file, self.id, ligand_name) f.write(dock_text) os.chdir(self.dock_dir) cmd = "dock6 -i {} -o {} -v".format(in_name, out_name) os.system(cmd) if ssbio.utils.is_non_zero_file(ranked_out): self.dock_flexible_outfile = out_name self.dock_flexible_conformers_result = conformers_out self.dock_flexible_scored_result = scored_out log.debug('{}: successful docking!'.format(self.dock_flexible_outfile)) else: log.error('{}+{}: empty DOCK6 ranked file, execution error (or ligand failed to dock)'.format(self.id,
op.basename(ligand_path)))
[docs] def auto_flexdock(self, binding_residues, radius, ligand_path=None, force_rerun=False): """Run DOCK6 on a PDB file, given its binding residues and a radius around them. Provide a path to a ligand to dock a ligand to it. If no ligand is provided, DOCK6 preparations will be run on that structure file. Args: binding_residues (str): Comma separated string of residues (eg: '144,170,199') radius (int, float): Radius around binding residues to dock to ligand_path (str): Path to ligand (mol2 format) to dock to protein force_rerun (bool): If method should be rerun even if output files exist """ log.debug('\n{}: running DOCK6...\n' '\tBinding residues: {}\n' '\tBinding residues radius: {}\n' '\tLigand to dock: {}\n'.format(self.id, binding_residues, radius, op.basename(ligand_path))) self.dockprep(force_rerun=force_rerun) self.protein_only_and_noH(force_rerun=force_rerun) self.dms_maker(force_rerun=force_rerun) self.sphgen(force_rerun=force_rerun) self.binding_site_mol2(residues=binding_residues, force_rerun=force_rerun) self.sphere_selector_using_residues(radius=radius, force_rerun=force_rerun) self.showbox(force_rerun=force_rerun) self.grid(force_rerun=force_rerun) if ligand_path:
self.do_dock6_flexible(ligand_path=ligand_path, force_rerun=force_rerun)
[docs]def parse_results_mol2(mol2_outpath): """Parse a DOCK6 mol2 output file, return a Pandas DataFrame of the results. Args: mol2_outpath (str): Path to mol2 output file Returns: DataFrame: Pandas DataFrame of the results """ docked_ligands = pd.DataFrame() lines = [line.strip() for line in open(mol2_outpath, 'r')] props = {} for i, line in enumerate(lines): if line.startswith('########## Name:'): ligand = line.strip().strip('##########').replace(' ', '').replace('\t', '').split(':')[1] line = lines[i + 1] props = {} props['Ligand'] = ligand if line.startswith('##########'): splitter = line.strip().strip('##########').replace(' ', '').replace('\t', '').split(':') props[splitter[0]] = float(splitter[1]) if line.startswith('@<TRIPOS>MOLECULE'): if props: docked_ligands = docked_ligands.append(props, ignore_index=True)
return docked_ligands # def do_dock6_rigid(self, ligand_path, force_rerun=False): # ligand_name = os.path.basename(args.ligand).split('.')[0] # in_name = op.join(self.dock_dir, "{}_{}_dock.in".format(self.id, ligand_name)) # out_name = op.join(self.dock_dir, "{}_{}_dock.out".format (self.id, ligand_name)) # # with open(in_name, "w") as f: # dock_text = """ligand_atom_file {} # limit_max_ligands no # skip_molecule no # read_mol_solvation no # calculate_rmsd no # use_database_filter no # orient_ligand yes # automated_matching yes # receptor_site_file {} # max_orientations 1000 # critical_points no # chemical_matching no # use_ligand_spheres no # use_internal_energy yes # internal_energy_rep_exp 12 # flexible_ligand no # bump_filter no # score_molecules yes # contact_score_primary no # contact_score_secondary no # grid_score_primary yes # grid_score_secondary no # grid_score_rep_rad_scale 1 # grid_score_vdw_scale 1 # grid_score_es_scale 1 # grid_score_grid_prefix {} # multigrid_score_secondary no # dock3.5_score_secondary no # continuous_score_secondary no # descriptor_score_secondary no # gbsa_zou_score_secondary no # gbsa_hawkins_score_secondary no # SASA_descriptor_score_secondary no # amber_score_secondary no # minimize_ligand yes # simplex_max_iterations 1000 # simplex_tors_premin_iterations 0 # simplex_max_cycles 1 # simplex_score_converge 0.1 # simplex_cycle_converge 1.0 # simplex_trans_step 1.0 # simplex_rot_step 0.1 # simplex_tors_step 10.0 # simplex_random_seed 0 # simplex_restraint_min yes # simplex_coefficient_restraint 10.0 # atom_model all # vdw_defn_file {} # flex_defn_file {} # flex_drive_file {} # ligand_outfile_prefix {}_{}_rigid # write_orientations yes # num_scored_conformers 20 # rank_ligands no # """.format(ligand_path, self.sphsel_path, self.grid_path.split('.')[0], # self.amb_file, self.flex1_file, self.flex2_file, self.id, ligand_name) # # f.write(dock_text) # # cmd = "dock6 -i {} -o {}".format(in_name, out_name) # os.system(cmd) # def do_dock6_amberscore(self, ligand_path, force_rerun=False): # """INCOMPLETE""" # ligand_name = os.path.basename(args.ligand).split('.')[0] # in_name = op.join(self.dock_dir, "{}_{}_amberscore.in".format(self.id, ligand_name)) # out_name = op.join(self.dock_dir, "{}_{}_amberscore.out".format(self.id, ligand_name)) # # with open(in_name, "w") as f: # dock_text = """ligand_atom_file {}.amber_score.mol2 # limit_max_ligands no # skip_molecule no # read_mol_solvation no # calculate_rmsd no # use_database_filter no # orient_ligand no # use_internal_energy no # flexible_ligand no # bump_filter no # score_molecules yes # contact_score_primary no # contact_score_secondary no # grid_score_primary no # grid_score_secondary no # multigrid_score_primary no # multigrid_score_secondary no # dock3.5_score_primary no # dock3.5_score_secondary no # continuous_score_primary no # continuous_score_secondary no # descriptor_score_primary no # descriptor_score_secondary no # gbsa_zou_score_primary no # gbsa_zou_score_secondary no # gbsa_hawkins_score_primary no # gbsa_hawkins_score_secondary no # SASA_descriptor_score_primary no # SASA_descriptor_score_secondary no # amber_score_primary yes # amber_score_secondary no # amber_score_receptor_file_prefix {} # amber_score_movable_region ligand # amber_score_minimization_rmsgrad 0.01 # amber_score_before_md_minimization_cycles 100 # amber_score_md_steps 3000 # amber_score_after_md_minimization_cycles 100 # amber_score_gb_model 5 # amber_score_nonbonded_cutoff 18.0 # amber_score_temperature 300.0 # amber_score_abort_on_unprepped_ligand yes # ligand_outfile_prefix output # write_orientations no # num_scored_conformers 1 # rank_ligands no # """.format() # # f.write(dock_text) # # cmd = "dock6 -i {} -o {} -v".format(in_name, out_name) # os.system(cmd) # if __name__ == '__main__': # # import glob # import argparse # import shlex # # # load inputs from command line # p = argparse.ArgumentParser( # description='Run the DOCK steps on a folder of structure files. To run in the background, execute using nohup: nohup dock.py $BASENAME $NUMFRAMES /path/to/structures/ /path/to/parameters/ --ligand /path/to/ligand.mol2 --cofactors $COFACTORS --residues $RESIDUES > /path/to/logs/$LOGNAME &') # p.add_argument('basename', help='base filename that you used to name your files') # p.add_argument('numframes', help='total number of frames from your trajectory', type=int) # p.add_argument('folder', help='path to folder with your structure files') # p.add_argument('params', help='path to folder with parameter files') # p.add_argument('--ligand', help='path to file of your ligand that you want to dock') # p.add_argument('--cofactors', # help='comma-separated string of cofactors that you want to keep while docking (e.g. SAM,DNC,WAT)') # p.add_argument('--residues', help='comma-separated string of the binding residues') # p.add_argument('--radius', help='radius around binding residues to dock to (default 9 A)', type=int, default=9) # p.add_argument('--redock', help='run DOCK again for the specified ligand, even if docking files exist', # default=False) # args = p.parse_args() # # # file paths for docking parameters # amb = os.path.join(args.params, 'vdw_AMBER_parm99.defn') # f1 = os.path.join(args.params, 'flex.defn') # f2 = os.path.join(args.params, 'flex_drive.tbl') # # print(args) # # loading current files # os.chdir(args.folder) # # pdbs = glob.glob('{}-*.pdb'.format(args.basename)) # current_files = os.listdir(os.getcwd()) # # # ligand name # if args.ligand: # ligandname = os.path.basename(args.ligand) # # # cofactors # if args.cofactors: # cofactors_list = shlex.split(args.cofactors) # else: # cofactors_list = [] # # print('***************PARAMETERS***************') # print('FULL LIST: {0}'.format(vars(args))) # if args.ligand: # print('LIGAND: {0}'.format(ligandname)) # if args.cofactors: # print('COFACTORS: {0}'.format(cofactors_list)) # if args.residues: # print('BINDING RESIDUES: {0}'.format(args.residues)) # print('RADIUS: {0}'.format(args.radius)) # # counter = 1 # for frame in range(1, args.numframes + 1): # # just a simple counter # print(str(counter) + '/' + str(args.numframes)) # counter += 1 # # # file_prefix = '{0}-{1:03d}'.format(args.basename, frame) # file_prefix = '{0}'.format(args.basename) # print(file_prefix) # # # DOCKPREP # # example: 3bwm-440_prep.mol2 # pdb = '{0}.pdb'.format(file_prefix) # prepped_check = '{}_prep.mol2'.format(file_prefix) # if prepped_check in current_files: # print('***DOCKPREP PREVIOUSLY RUN***') # prepped_file = prepped_check # else: # print('RUNNING: DOCKPREP') # prepped_file = dockprep(pdb, file_prefix) # # # ISOLATE RECEPTOR # # example: 3bwm-440_receptor.mol2, 3bwm-440_receptor_noH.pdb # receptor_check = '{}_receptor.mol2'.format(file_prefix) # receptor_noH_check = '{}_receptor_noH.pdb'.format(file_prefix) # if receptor_check in current_files and receptor_noH_check in current_files: # print('***RECEPTOR FILES PREVIOUSLY GENERATED***') # receptor, receptor_noH = receptor_check, receptor_noH_check # else: # print('RUNNING: ISOLATE RECEPTOR') # receptor, receptor_noH = protein_only_and_noH(prepped_file, cofactors_list, file_prefix) # # # DMS # # example: 3bwm-440_receptor.dms # dms_check = '{}_receptor.dms'.format(file_prefix) # if dms_check in current_files: # print('***SURFACE PREVIOUSLY GENERATED***') # dms = dms_check # else: # print('RUNNING: DMS') # dms = dms_maker(receptor_noH, file_prefix) # # # SPHGEN # # example: 3bwm-440_receptor.sph # sph_check = '{}_receptor.sph'.format(file_prefix) # if sph_check in current_files: # print('***SPHERES PREVIOUSLY GENERATED***') # sph = sph_check # else: # print('RUNNING: SPHGEN') # sph = sphgen(dms, file_prefix) # # # SPHERE_SELECTOR # # first choose binding site and save it as separate .mol2 # # example: 3BWY-418_binding_residues.mol2 # binding_site_mol2_file_check = '{}_binding_residues.mol2'.format(file_prefix) # if binding_site_mol2_file_check in current_files: # print('***BINDING SITE RESIDUES ALREADY DEFINED***') # binding_site_mol2_file = binding_site_mol2_file_check # else: # print('RUNNING: BINDING SITE MOL2') # binding_site_mol2_file = binding_site_mol2(receptor_noH, args.residues, file_prefix) # # # then select the spheres based on these binding residues # # example: 3bwm-440_selected_spheres_using_binding_residues.sph # sel_sph_check = '{}_selected_spheres_using_binding_residues.sph'.format(file_prefix) # if sel_sph_check in current_files: # print('***SPHERES ALREADY SELECTED***') # sel_sph = sel_sph_check # else: # print('RUNNING: SPHERE_SELECTOR') # sel_sph = sphere_selector_using_residues(sph, binding_site_mol2_file, args.radius, file_prefix) # # # SHOWBOX # # example: 3bwm-440_box.pdb # box_check = '{}_box.pdb'.format(file_prefix) # if box_check in current_files: # print('***BOX PREVIOUSLY MADE***') # box = box_check # else: # print('RUNNING: SHOWBOX') # box = showbox(sel_sph, file_prefix) # # # GRID # # example: 3bwm-440_grid.out # gr_check = '{}_grid.out'.format(file_prefix) # if gr_check in current_files: # print('***GRID PREVIOUSLY CALCULATED***') # gr = gr_check # else: # print('RUNNING: GRID') # gr = grid(receptor, box, amb, file_prefix) # # # DOCK # if args.ligand: # dock6_flexible_check = '{}_{}_flexible_scored.mol2'.format((file_prefix, ligandname.split('.')[0])) # if dock6_flexible_check in current_files and not args.redock: # print('***DOCK PREVIOUSLY RUN***') # else: # print('RUNNING: DOCK') # do_dock6_flexible(args.ligand, sel_sph, gr, amb, f1, f2, file_prefix) # # print('***DOCKING COMPLETE***')