diff --git a/bin/post/lgd_generate_conformations.py b/bin/post/lgd_generate_conformations.py index 0d11cfb..bf8d79b 100755 --- a/bin/post/lgd_generate_conformations.py +++ b/bin/post/lgd_generate_conformations.py @@ -13,6 +13,7 @@ from lightdock.structure.complex import Complex from lightdock.mathutil.cython.quaternion import Quaternion from lightdock.structure.nm import read_nmodes +from lightdock.prep.simulation import get_setup_from_file log = LoggingManager.get_logger('generate_conformations') @@ -59,7 +60,7 @@ def get_lightdock_structures(input_file): return file_names -def parse_output_file(lightdock_output): +def parse_output_file(lightdock_output, num_anm_rec, num_anm_lig): translations = [] rotations = [] receptor_ids = [] @@ -79,9 +80,9 @@ def parse_output_file(lightdock_output): coord = line[1:last].split(',') translations.append([float(coord[0]), float(coord[1]), float(coord[2])]) rotations.append(Quaternion(float(coord[3]), float(coord[4]), float(coord[5]), float(coord[6]))) - if len(coord) == (7 + DEFAULT_NMODES_REC + DEFAULT_NMODES_LIG): - rec_extents.append(np.array([float(x) for x in coord[7:7+DEFAULT_NMODES_REC]])) - lig_extents.append(np.array([float(x) for x in coord[-DEFAULT_NMODES_LIG:]])) + if len(coord) > 7: + rec_extents.append(np.array([float(x) for x in coord[7:7+num_anm_rec]])) + lig_extents.append(np.array([float(x) for x in coord[-num_anm_lig:]])) raw_data = line[last+1:].split() receptor_id = int(raw_data[0]) ligand_id = int(raw_data[1]) @@ -105,9 +106,21 @@ def parse_output_file(lightdock_output): type=valid_file, metavar="lightdock_output") # Number of glowworms parser.add_argument("glowworms", help="number of glowworms", type=valid_integer_number) + # Optional, setup file + parser.add_argument("--setup", "-setup", "-s", help="Simulation setup file", + dest="setup_file", metavar="setup_file", type=valid_file, default=None) args = parser.parse_args() + # Load setup configuration if provided + setup = get_setup_from_file(args.setup_file) if args.setup_file else None + + num_anm_rec = DEFAULT_NMODES_REC + num_anm_lig = DEFAULT_NMODES_LIG + if setup and setup['use_anm']: + num_anm_rec = setup['anm_rec'] + num_anm_lig = setup['anm_lig'] + # Receptor structures = [] for structure in get_lightdock_structures(args.receptor_structures): @@ -128,7 +141,7 @@ def parse_output_file(lightdock_output): # Output file translations, rotations, receptor_ids, ligand_ids, \ - rec_extents, lig_extents = parse_output_file(args.lightdock_output) + rec_extents, lig_extents = parse_output_file(args.lightdock_output, num_anm_rec, num_anm_lig) found_conformations = len(translations) num_conformations = args.glowworms @@ -157,22 +170,30 @@ def parse_output_file(lightdock_output): # Use normal modes if provided: if len(rec_extents): try: - for nm in range(DEFAULT_NMODES_REC): + for nm in range(num_anm_rec): receptor_pose.coordinates += nmodes_rec[nm] * rec_extents[i][nm] - except ValueError, e: + except ValueError: log.error("Problem found on calculating ANM for receptor:") log.error("Number of atom coordinates is: %s" % str(receptor_pose.coordinates.shape)) log.error("Number of ANM is: %s" % str(nmodes_rec.shape)) raise SystemExit + except IndexError: + log.error("Problem found on calculating ANM for receptor:") + log.error("If you have used anm_rec different than default, please use --setup") + raise SystemExit if len(lig_extents): try: - for nm in range(DEFAULT_NMODES_LIG): + for nm in range(num_anm_lig): ligand_pose.coordinates += nmodes_lig[nm] * lig_extents[i][nm] - except ValueError, e: + except ValueError: log.error("Problem found on calculating ANM for ligand:") log.error("Number of atom coordinates is: %s" % str(receptor_pose.coordinates.shape)) log.error("Number of ANM is: %s" % str(nmodes_rec.shape)) raise SystemExit + except IndexError: + log.error("Problem found on calculating ANM for ligand:") + log.error("If you have used anm_lig different than default, please use --setup") + raise SystemExit # We rotate first, ligand it's at initial position ligand_pose.rotate(rotations[i]) diff --git a/bin/simulation/docking_mpi.py b/bin/simulation/docking_mpi.py index f680055..e20a0b7 100644 --- a/bin/simulation/docking_mpi.py +++ b/bin/simulation/docking_mpi.py @@ -38,6 +38,8 @@ def set_gso(number_of_glowworms, adapters, scoring_functions, initial_positions, else: gso_parameters = GSOParameters() builder = LightdockGSOBuilder() + if not use_anm: + anm_rec = anm_lig = 0 gso = builder.create_from_file(number_of_glowworms, random_number_generator, gso_parameters, adapters, scoring_functions, bounding_box, initial_positions, step_translation, step_rotation, nmodes_step, local_minimization) @@ -112,8 +114,16 @@ def run_simulation(parser): ligand.move_to_origin() if args.use_anm: - receptor.n_modes = read_nmodes("%s%s" % (DEFAULT_REC_NM_FILE, NUMPY_FILE_SAVE_EXTENSION) ) - ligand.n_modes = read_nmodes("%s%s" % (DEFAULT_LIG_NM_FILE, NUMPY_FILE_SAVE_EXTENSION) ) + try: + receptor.n_modes = read_nmodes("%s%s" % (DEFAULT_REC_NM_FILE, NUMPY_FILE_SAVE_EXTENSION) ) + except: + log.warning("No ANM found for receptor molecule") + receptor.n_modes = None + try: + ligand.n_modes = read_nmodes("%s%s" % (DEFAULT_LIG_NM_FILE, NUMPY_FILE_SAVE_EXTENSION) ) + except: + log.warning("No ANM found for ligand molecule") + ligand.n_modes = None starting_points_files = load_starting_positions(args.swarms, args.glowworms, args.use_anm, args.anm_rec, args.anm_lig) diff --git a/bin/simulation/docking_multiprocessing.py b/bin/simulation/docking_multiprocessing.py index 258c1b1..8008535 100644 --- a/bin/simulation/docking_multiprocessing.py +++ b/bin/simulation/docking_multiprocessing.py @@ -37,9 +37,12 @@ def set_gso(number_of_glowworms, adapters, scoring_functions, initial_positions, else: gso_parameters = GSOParameters() builder = LightdockGSOBuilder() + if not use_anm: + anm_rec = anm_lig = 0 gso = builder.create_from_file(number_of_glowworms, random_number_generator, gso_parameters, adapters, scoring_functions, bounding_box, initial_positions, - step_translation, step_rotation, nmodes_step, local_minimization) + step_translation, step_rotation, nmodes_step, local_minimization, + anm_rec, anm_lig) return gso @@ -121,8 +124,16 @@ def run_simulation(parser): ligand.move_to_origin() if args.use_anm: - receptor.n_modes = read_nmodes("%s%s" % (DEFAULT_REC_NM_FILE, NUMPY_FILE_SAVE_EXTENSION) ) - ligand.n_modes = read_nmodes("%s%s" % (DEFAULT_LIG_NM_FILE, NUMPY_FILE_SAVE_EXTENSION) ) + try: + receptor.n_modes = read_nmodes("%s%s" % (DEFAULT_REC_NM_FILE, NUMPY_FILE_SAVE_EXTENSION) ) + except: + log.warning("No ANM found for receptor molecule") + receptor.n_modes = None + try: + ligand.n_modes = read_nmodes("%s%s" % (DEFAULT_LIG_NM_FILE, NUMPY_FILE_SAVE_EXTENSION) ) + except: + log.warning("No ANM found for ligand molecule") + ligand.n_modes = None starting_points_files = load_starting_positions(args.swarms, args.glowworms, args.use_anm, args.anm_rec, args.anm_lig) diff --git a/bin/simulation/lightdock_setup.py b/bin/simulation/lightdock_setup.py index a96f7f6..31d37ec 100755 --- a/bin/simulation/lightdock_setup.py +++ b/bin/simulation/lightdock_setup.py @@ -61,8 +61,12 @@ # Calculate and save ANM if required if args.use_anm: - calculate_anm(receptor, args.anm_rec, DEFAULT_REC_NM_FILE) - calculate_anm(ligand, args.anm_lig, DEFAULT_LIG_NM_FILE) + if args.anm_rec > 0: + log.info("Calculating ANM for receptor molecule...") + calculate_anm(receptor, args.anm_rec, DEFAULT_REC_NM_FILE) + if args.anm_lig > 0: + log.info("Calculating ANM for ligand molecule...") + calculate_anm(ligand, args.anm_lig, DEFAULT_LIG_NM_FILE) # Parse restraints if any: receptor_restraints = ligand_restraints = None diff --git a/lightdock/gso/algorithm.py b/lightdock/gso/algorithm.py index ade81f8..8370c8c 100644 --- a/lightdock/gso/algorithm.py +++ b/lightdock/gso/algorithm.py @@ -107,13 +107,14 @@ class LightdockGSOBuilder(object): """Creates a GSO simulation object for the LightDock framework""" def create_from_file(self, number_of_glowworms, random_number_generator, gso_parameters, adapters, scoring_functions, bounding_box, initial_population_file, - step_translation, step_rotation, step_nmodes, local_minimization): + step_translation, step_rotation, step_nmodes, local_minimization, + anm_rec, anm_lig): """Creates a new GSO instance of the algorithm reading the initial position of the glowworms agents from initial_population_file and using the scoring function adapter. """ self._initializer = LightdockFromFileInitializer(adapters, scoring_functions, number_of_glowworms, gso_parameters, bounding_box.dimension, initial_population_file, step_translation, step_rotation, - random_number_generator, step_nmodes) + random_number_generator, step_nmodes, anm_rec, anm_lig) return GSO(self._initializer.generate_glowworms(), gso_parameters, random_number_generator, local_minimization=local_minimization) diff --git a/lightdock/gso/initializer.py b/lightdock/gso/initializer.py index 0c778a3..e9c2dfa 100644 --- a/lightdock/gso/initializer.py +++ b/lightdock/gso/initializer.py @@ -88,7 +88,7 @@ class LightdockFromFileInitializer(Initializer): """This initializer takes into account the complex provided by the adapter""" def __init__(self, adapters, scoring_functions, number_of_glowworms, gso_parameters, dimensions, initial_population_file, step_translation, step_rotation, - random_number_generator, step_nmodes): + random_number_generator, step_nmodes, anm_rec, anm_lig): super(LightdockFromFileInitializer, self).__init__(scoring_functions, number_of_glowworms, gso_parameters) self.dimensions = dimensions self.initial_population_file = initial_population_file @@ -98,6 +98,8 @@ def __init__(self, adapters, scoring_functions, number_of_glowworms, gso_paramet self.step_nmodes = step_nmodes # Patch to not mess with old simulations self.random_number_generator = MTGenerator(random_number_generator.seed) + self.anm_rec = anm_rec + self.anm_lig = anm_lig def generate_landscape_positions(self): """Generates a list of landscape positions that have been read @@ -122,6 +124,5 @@ def generate_landscape_positions(self): adapter.receptor_model, adapter.ligand_model, receptor_index, ligand_index, self.step_translation, self.step_rotation, - self.step_nmodes)) - + self.step_nmodes, self.anm_rec, self.anm_lig)) return positions diff --git a/lightdock/gso/searchspace/landscape.py b/lightdock/gso/searchspace/landscape.py index 0d345b0..b592d69 100644 --- a/lightdock/gso/searchspace/landscape.py +++ b/lightdock/gso/searchspace/landscape.py @@ -97,17 +97,10 @@ class DockingLandscapePosition(LandscapePosition): """ def __init__(self, scoring_function, coordinates, receptor, ligand, receptor_id=0, ligand_id=0, step_translation=DEFAULT_TRANSLATION_STEP, step_rotation=DEFAULT_ROTATION_STEP, - step_nmodes=DEFAULT_NMODES_STEP): + step_nmodes=0, num_rec_nmodes=0, num_lig_nmodes=0): self.objective_function = scoring_function self.translation = np.array(coordinates[:3]) self.rotation = Quaternion(coordinates[3], coordinates[4], coordinates[5], coordinates[6]) - # Copy ANM information if required - if len(coordinates) > 7 and len(coordinates) == (7 + DEFAULT_NMODES_REC + DEFAULT_NMODES_LIG): - self.rec_extent = np.array(coordinates[7:7+DEFAULT_NMODES_REC]) - self.lig_extent = np.array(coordinates[-DEFAULT_NMODES_LIG:]) - else: - self.rec_extent = np.array([]) - self.lig_extent = np.array([]) self.receptor = receptor self.ligand = ligand self.receptor_id = receptor_id @@ -115,6 +108,11 @@ def __init__(self, scoring_function, coordinates, receptor, ligand, receptor_id= self.step_translation = step_translation self.step_rotation = step_rotation self.step_nmodes = step_nmodes + self.num_rec_nmodes = num_rec_nmodes + self.num_lig_nmodes = num_lig_nmodes + # Copy ANM information if required + self.rec_extent = np.array(coordinates[7:7+self.num_rec_nmodes]) if self.num_rec_nmodes > 0 else np.array([]) + self.lig_extent = np.array(coordinates[-self.num_lig_nmodes:]) if self.num_lig_nmodes > 0 else np.array([]) # This part is important, each position needs to retain its own pose coordinates self.receptor_pose = self.receptor.coordinates[self.receptor_id].clone() self.ligand_pose = self.ligand.coordinates[self.ligand_id].clone() @@ -129,7 +127,8 @@ def clone(self): return DockingLandscapePosition(self.objective_function, coordinates, self.receptor, self.ligand, self.receptor_id, self.ligand_id, - self.step_translation, self.step_rotation, self.step_nmodes) + self.step_translation, self.step_rotation, self.step_nmodes, + self.num_rec_nmodes, self.num_lig_nmodes) def evaluate_objective_function(self, receptor_structure_id=None, ligand_structure_id=None): """Evaluates the objective function at the given coordinates""" @@ -147,13 +146,11 @@ def evaluate_objective_function(self, receptor_structure_id=None, ligand_structu self.ligand_reference_points = self.ligand.reference_points.clone() # Use normal modes if provided: - num_rec_nmodes = len(self.rec_extent) - if num_rec_nmodes: - for i in range(num_rec_nmodes): + if self.num_rec_nmodes > 0: + for i in range(self.num_rec_nmodes): self.receptor_pose.coordinates += self.receptor.n_modes[i] * self.rec_extent[i] - num_lig_nmodes = len(self.lig_extent) - if num_lig_nmodes: - for i in range(num_lig_nmodes): + if self.num_lig_nmodes > 0: + for i in range(self.num_lig_nmodes): self.ligand_pose.coordinates += self.ligand.n_modes[i] * self.lig_extent[i] # We rotate first, ligand it's at initial position @@ -198,14 +195,14 @@ def move(self, other): # Rotation (Quaternion SLERP) self.rotation = self.rotation.slerp(other.rotation, self.step_rotation) # NModes - if len(self.rec_extent): + if self.num_rec_nmodes > 0: delta_x = other.rec_extent - self.rec_extent n = np.linalg.norm(delta_x) # Only move if required if not np.allclose([0.0], [n]): delta_x *= (self.step_nmodes / n) self.rec_extent += delta_x - if len(self.lig_extent): + if self.num_lig_nmodes > 0: delta_x = other.lig_extent - self.lig_extent n = np.linalg.norm(delta_x) # Only move if required @@ -238,8 +235,8 @@ def update_landscape_position(self, optimized_vector): """Updates the current pose""" self.translation = optimized_vector[:3] self.rotation = Quaternion(optimized_vector[3], optimized_vector[4], optimized_vector[5], optimized_vector[6]) - self.rec_extent = optimized_vector[7:7+DEFAULT_NMODES_REC] - self.lig_extent = optimized_vector[-DEFAULT_NMODES_LIG:] + self.rec_extent = optimized_vector[7:7+self.num_rec_nmodes] if self.num_rec_nmodes > 0 else np.array([]) + self.lig_extent = optimized_vector[-self.num_lig_nmodes:] if self.num_lig_nmodes > 0 else np.array([]) def minimize(self): """Returns the new scoring after minimizing this landscape position using a local non-grandient diff --git a/lightdock/test/gso/test_initializer.py b/lightdock/test/gso/test_initializer.py index 29b0a50..29cd36b 100644 --- a/lightdock/test/gso/test_initializer.py +++ b/lightdock/test/gso/test_initializer.py @@ -109,7 +109,7 @@ def test_create_swarm(self): initializer = LightdockFromFileInitializer([self.adapter], [self.scoring_function], number_of_glowworms, gso_parameters, 7, self.golden_data_path+'initial_positions_1PPE.txt', - 0.5, 0.5, random_number_generator, 0.5) + 0.5, 0.5, random_number_generator, 0.5, 10, 10) swarm = initializer.generate_glowworms() assert number_of_glowworms == swarm.get_size() @@ -123,7 +123,7 @@ def test_generate_landscape_positions_without_coordinates(self): initializer = LightdockFromFileInitializer(self.adapter, self.scoring_function, number_of_glowworms, gso_parameters, 7, self.golden_data_path+'initial_positions_empty.txt', - 0.5, 0.5, random_number_generator, 0.5) + 0.5, 0.5, random_number_generator, 0.5, 10, 10) swarm = initializer.generate_glowworms() assert swarm.get_size() > 0 @@ -137,7 +137,7 @@ def test_generate_landscape_positions_num_glowworms_different(self): initializer = LightdockFromFileInitializer(self.adapter, self.scoring_function, number_of_glowworms, gso_parameters, 7, self.golden_data_path+'initial_positions_1PPE.txt', - 0.5, 0.5, random_number_generator, 0.5) + 0.5, 0.5, random_number_generator, 0.5, 10, 10) swarm = initializer.generate_glowworms() assert swarm.get_size() > 0 diff --git a/lightdock/test/gso/test_lightdockbuilder.py b/lightdock/test/gso/test_lightdockbuilder.py index 87816b2..bc9007f 100644 --- a/lightdock/test/gso/test_lightdockbuilder.py +++ b/lightdock/test/gso/test_lightdockbuilder.py @@ -53,7 +53,7 @@ def test_LightDockGSOBuilder_using_FromFileInitializer(self): gso = builder.create_from_file(number_of_glowworms, self.random_number_generator, self.gso_parameters, [adapter], [scoring_function], self.bounding_box, self.golden_data_path+'initial_positions_1PPE.txt', - 0.5, 0.5, 0.5, local_minimization=False) + 0.5, 0.5, 0.5, False, 10, 10) assert 5 == gso.swarm.get_size() diff --git a/lightdock/util/parser.py b/lightdock/util/parser.py index 001e6b1..fd26f4f 100644 --- a/lightdock/util/parser.py +++ b/lightdock/util/parser.py @@ -49,6 +49,16 @@ def valid_integer_number(int_value): if int_value <= 0: raise argparse.ArgumentTypeError("%s is an invalid value" % int_value) return int_value + + @staticmethod + def valid_natural_number(int_value): + try: + int_value = int(int_value) + except: + raise argparse.ArgumentTypeError("%s is an invalid value" % int_value) + if int_value < 0: + raise argparse.ArgumentTypeError("%s is an invalid value" % int_value) + return int_value @staticmethod def valid_float_number(float_value): @@ -92,10 +102,10 @@ def __init__(self): parser.add_argument("--seed_anm", help="Random seed used in ANM intial extent", dest="anm_seed", type=int, default=STARTING_NM_SEED) parser.add_argument("-anm_rec", "--anm_rec", help="Number of ANM modes for receptor", - type=SetupCommandLineParser.valid_integer_number, + type=SetupCommandLineParser.valid_natural_number, dest="anm_rec", default=DEFAULT_NMODES_REC) parser.add_argument("-anm_lig", "--anm_lig", help="Number of ANM modes for ligand", - type=SetupCommandLineParser.valid_integer_number, + type=SetupCommandLineParser.valid_natural_number, dest="anm_lig", default=DEFAULT_NMODES_LIG) # Restraints file parser.add_argument("-rst", "--rst", help="Restraints file", diff --git a/lightdock/version.py b/lightdock/version.py index 17399c8..45d0047 100644 --- a/lightdock/version.py +++ b/lightdock/version.py @@ -1,3 +1,3 @@ """Framework version""" -CURRENT_VERSION = "0.5.6" +CURRENT_VERSION = "0.6.0"