From 27c55e352161c7b05f2b592296441b3e0940a601 Mon Sep 17 00:00:00 2001 From: karolamik13 Date: Sat, 14 Oct 2023 22:32:33 +0200 Subject: [PATCH 01/11] calcLigandInteractions() - tempfile & improvements --- prody/proteins/interactions.py | 39 +++++++++++++++++----------------- 1 file changed, 19 insertions(+), 20 deletions(-) diff --git a/prody/proteins/interactions.py b/prody/proteins/interactions.py index 95d1d7067..13bd84d42 100644 --- a/prody/proteins/interactions.py +++ b/prody/proteins/interactions.py @@ -2021,9 +2021,9 @@ def calcLigandInteractions(atoms, **kwargs): :arg atoms: an Atomic object from which residues are selected :type atoms: :class:`.Atomic` - :arg select: a selection string for residues of interest + :arg sele: a selection string for residues of interest default is 'all not (water or protein or ion)' - :type select: str + :type sele: str :arg ignore_ligs: List of ligands which will be excluded from the analysis. :type ignore_ligs: list @@ -2056,34 +2056,33 @@ def calcLigandInteractions(atoms, **kwargs): 'with `getCoords` method') try: from plip.structure.preparation import PDBComplex - - pdb_name = atoms.getTitle()+'_sele.pdb' - LOGGER.info("Writing PDB file with selection in the local directory.") - writePDB(pdb_name, atoms) + import tempfile + + with tempfile.NamedTemporaryFile(delete=False, suffix=".pdb") as temp_pdb_file: + writePDB(temp_pdb_file.name, atoms) + temp_pdb_file_name = temp_pdb_file.name + + #pdb_name = atoms.getTitle()+'_sele.pdb' + #LOGGER.info("Writing PDB file with selection in the local directory.") + #writePDB(pdb_name, atoms) try: if atoms.hydrogen == None or atoms.hydrogen.numAtoms() < 30: # if there is no hydrogens in PDB structure - addMissingAtoms(pdb_name) - pdb_name = pdb_name[:-4]+'_addH.pdb' - atoms = parsePDB(pdb_name) - LOGGER.info("Lack of hydrogens in the structure. Hydrogens have been added.") + LOGGER.info("Lack of hydrogens in the structure. Use addMissingAtoms to add them.") except: - raise ValueError("Missing atoms from the side chains of the structure. Use addMissingAtoms.") - + pass + Ligands = [] # Ligands can be more than one my_mol = PDBComplex() - my_mol.load_pdb(pdb_name) # Load the PDB file into PLIP class + my_mol.load_pdb(temp_pdb_file_name) # Load the PDB file into PLIP class - if 'select' in kwargs: - select = kwargs['select'] - LOGGER.info('Selection will be replaced.') + if 'sele' in kwargs: + sele = kwargs['sele'] else: select='all not (water or protein or ion)' - LOGGER.info('Default selection will be used.') if 'ignore_ligs' in kwargs: ignore_ligs = kwargs['ignore_ligs'] - LOGGER.info('Ignoring list of ligands is uploaded.') else: ignore_ligs=['NAG','BMA','MAN'] LOGGER.info('Three molecules will be ignored from analysis: NAG, BMA and MAN.') @@ -2109,8 +2108,8 @@ def calcLigandInteractions(atoms, **kwargs): return Ligands, analyzedLigand - except: - LOGGER.info("Install Openbabel and PLIP.") + except ImportError: + raise ImportError("Install Openbabel and PLIP.") def listLigandInteractions(PLIP_output): From 67155244eb0b3cfe75a37c5cad801e54b97130b7 Mon Sep 17 00:00:00 2001 From: karolamik13 Date: Sun, 15 Oct 2023 14:56:48 +0200 Subject: [PATCH 02/11] small changes in calcLigandInteractions --- prody/proteins/interactions.py | 7 +------ 1 file changed, 1 insertion(+), 6 deletions(-) diff --git a/prody/proteins/interactions.py b/prody/proteins/interactions.py index 13bd84d42..449c8b823 100644 --- a/prody/proteins/interactions.py +++ b/prody/proteins/interactions.py @@ -2062,10 +2062,6 @@ def calcLigandInteractions(atoms, **kwargs): writePDB(temp_pdb_file.name, atoms) temp_pdb_file_name = temp_pdb_file.name - #pdb_name = atoms.getTitle()+'_sele.pdb' - #LOGGER.info("Writing PDB file with selection in the local directory.") - #writePDB(pdb_name, atoms) - try: if atoms.hydrogen == None or atoms.hydrogen.numAtoms() < 30: # if there is no hydrogens in PDB structure LOGGER.info("Lack of hydrogens in the structure. Use addMissingAtoms to add them.") @@ -2084,8 +2080,7 @@ def calcLigandInteractions(atoms, **kwargs): if 'ignore_ligs' in kwargs: ignore_ligs = kwargs['ignore_ligs'] else: - ignore_ligs=['NAG','BMA','MAN'] - LOGGER.info('Three molecules will be ignored from analysis: NAG, BMA and MAN.') + ignore_ligs=['MAN'] select = select+' and not (resname '+' '.join(ignore_ligs)+')' ligand_select = atoms.select(select) From 2ec7c754362e5c6b7b3a0f713cefb5268fa1ff90 Mon Sep 17 00:00:00 2001 From: karolamik13 Date: Sun, 15 Oct 2023 15:36:15 +0200 Subject: [PATCH 03/11] reorganization and improvements in protein-ligand interactions --- prody/proteins/interactions.py | 112 +++++++++++++++++---------------- 1 file changed, 57 insertions(+), 55 deletions(-) diff --git a/prody/proteins/interactions.py b/prody/proteins/interactions.py index 449c8b823..0c5f3f032 100644 --- a/prody/proteins/interactions.py +++ b/prody/proteins/interactions.py @@ -41,7 +41,7 @@ 'calcPiCationTrajectory', 'calcHydrophobicTrajectory', 'calcDisulfideBondsTrajectory', 'calcProteinInteractions', 'calcStatisticsInteractions', 'calcDistribution', 'calcSASA', 'calcVolume','compareInteractions', 'showInteractionsGraph', - 'calcLigandInteractions', 'listLigandInteractions', + 'calcLigandInteractions', 'showLigandInteractions', 'showProteinInteractions_VMD', 'showLigandInteraction_VMD', 'calcHydrogenBondsTrajectory', 'calcHydrophobicOverlapingAreas', 'Interactions', 'InteractionsTrajectory'] @@ -2013,6 +2013,61 @@ def calcDistribution(interactions, residue1, residue2=None, **kwargs): for i in additional_residues: LOGGER.info(i) + +def showLigandInteractions(PLIP_output): + """Create a list of interactions from PLIP output created using calcLigandInteractions(). + Results can be displayed in VMD. + + :arg PLIP_output: Results from PLIP for protein-ligand interactions. + :type PLIP_output: PLIP object obtained from calcLigandInteractions() + + Note that five types of interactions are considered: hydrogen bonds, salt bridges, + pi-stacking, cation-pi, hydrophobic and water bridges.""" + + Inter_list_all = [] + for i in PLIP_output.all_itypes: + param_inter = [method for method in dir(i) if method.startswith('_') is False] + + if str(type(i)).split('.')[-1].strip("'>") == 'hbond': + Inter_list = ['hbond',i.restype+str(i.resnr), i[0].type+'_'+str(i.d_orig_idx), i.reschain, + i.restype+str(i.resnr_l), i[2].type+'_'+str(i.a_orig_idx), i.reschain_l, + i.distance_ad, i.angle, i[0].coords, i[2].coords] + + if str(type(i)).split('.')[-1].strip("'>") == 'saltbridge': + Inter_list = ['saltbridge',i.restype+str(i.resnr), '_'.join(map(str,i.negative.atoms_orig_idx)), i.reschain, + i.restype+str(i.resnr_l), '_'.join(map(str,i.positive.atoms_orig_idx)), i.reschain_l, + i.distance, None, i.negative.center, i.positive.center] + + if str(type(i)).split('.')[-1].strip("'>") == 'pistack': + Inter_list = ['pistack',i.restype+str(i.resnr), '_'.join(map(str,i[0].atoms_orig_idx)), i.reschain, + i.restype+str(i.resnr_l), '_'.join(map(str,i[1].atoms_orig_idx)), i.reschain_l, + i.distance, i.angle, i[0].center, i[1].center] + + if str(type(i)).split('.')[-1].strip("'>") == 'pication': + Inter_list = ['pication',i.restype+str(i.resnr), '_'.join(map(str,i[0].atoms_orig_idx)), i.reschain, + i.restype+str(i.resnr_l), '_'.join(map(str,i[1].atoms_orig_idx)), i.reschain_l, + i.distance, None, i[0].center, i[1].center] + + if str(type(i)).split('.')[-1].strip("'>") == 'hydroph_interaction': + Inter_list = ['hydroph_interaction',i.restype+str(i.resnr), i[0].type+'_'+str(i[0].idx), i.reschain, + i.restype+str(i.resnr_l), i[2].type+'_'+str(i[2].idx), i.reschain_l, + i.distance, None, i[0].coords, i[2].coords] + + if str(type(i)).split('.')[-1].strip("'>") == 'waterbridge': + water = i.water + Inter_list = ['waterbridge',i.restype+str(i.resnr), i[0].type+'_'+str(i[0].idx), i.reschain, + i.restype+str(i.resnr_l), i[3].type+'_'+str(i[3].idx), i.reschain_l, + [i.distance_aw, i.distance_dw], [i.d_angle, i.w_angle], i[0].coords, i[3].coords, + i.water.coords, i[7].residue.name+'_'+str(i[7].residue.idx)] + + Inter_list_all.append(Inter_list) + + for nr_k,k in enumerate(Inter_list_all): + LOGGER.info("%3i%22s%10s%26s%4s <---> %8s%12s%4s%6.1f" % (nr_k+1,k[0],k[1],k[2],k[3],k[4],k[5],k[6],k[7])) + + return Inter_list_all + + def calcLigandInteractions(atoms, **kwargs): """Provide ligand interactions with other elements of the system including protein, water and ions. Results are computed by PLIP [SS15]_ which should be installed. @@ -2098,6 +2153,7 @@ def calcLigandInteractions(atoms, **kwargs): my_mol.analyze() my_interactions = my_mol.interaction_sets[my_bsid] # Contains all interaction data Ligands.append(my_interactions) + showLigandInteractions(my_interactions) except: LOGGER.info(my_bsid+" not analyzed") @@ -2107,60 +2163,6 @@ def calcLigandInteractions(atoms, **kwargs): raise ImportError("Install Openbabel and PLIP.") -def listLigandInteractions(PLIP_output): - """Create a list of interactions from PLIP output created using calcLigandInteractions(). - Results can be displayed in VMD. - - :arg PLIP_output: Results from PLIP for protein-ligand interactions. - :type PLIP_output: PLIP object obtained from calcLigandInteractions() - - Note that five types of interactions are considered: hydrogen bonds, salt bridges, - pi-stacking, cation-pi, hydrophobic and water bridges.""" - - Inter_list_all = [] - for i in PLIP_output.all_itypes: - param_inter = [method for method in dir(i) if method.startswith('_') is False] - - if str(type(i)).split('.')[-1].strip("'>") == 'hbond': - Inter_list = ['hbond',i.restype+str(i.resnr), i[0].type+'_'+str(i.d_orig_idx), i.reschain, - i.restype+str(i.resnr_l), i[2].type+'_'+str(i.a_orig_idx), i.reschain_l, - i.distance_ad, i.angle, i[0].coords, i[2].coords] - - if str(type(i)).split('.')[-1].strip("'>") == 'saltbridge': - Inter_list = ['saltbridge',i.restype+str(i.resnr), '_'.join(map(str,i.negative.atoms_orig_idx)), i.reschain, - i.restype+str(i.resnr_l), '_'.join(map(str,i.positive.atoms_orig_idx)), i.reschain_l, - i.distance, None, i.negative.center, i.positive.center] - - if str(type(i)).split('.')[-1].strip("'>") == 'pistack': - Inter_list = ['pistack',i.restype+str(i.resnr), '_'.join(map(str,i[0].atoms_orig_idx)), i.reschain, - i.restype+str(i.resnr_l), '_'.join(map(str,i[1].atoms_orig_idx)), i.reschain_l, - i.distance, i.angle, i[0].center, i[1].center] - - if str(type(i)).split('.')[-1].strip("'>") == 'pication': - Inter_list = ['pication',i.restype+str(i.resnr), '_'.join(map(str,i[0].atoms_orig_idx)), i.reschain, - i.restype+str(i.resnr_l), '_'.join(map(str,i[1].atoms_orig_idx)), i.reschain_l, - i.distance, None, i[0].center, i[1].center] - - if str(type(i)).split('.')[-1].strip("'>") == 'hydroph_interaction': - Inter_list = ['hydroph_interaction',i.restype+str(i.resnr), i[0].type+'_'+str(i[0].idx), i.reschain, - i.restype+str(i.resnr_l), i[2].type+'_'+str(i[2].idx), i.reschain_l, - i.distance, None, i[0].coords, i[2].coords] - - if str(type(i)).split('.')[-1].strip("'>") == 'waterbridge': - water = i.water - Inter_list = ['waterbridge',i.restype+str(i.resnr), i[0].type+'_'+str(i[0].idx), i.reschain, - i.restype+str(i.resnr_l), i[3].type+'_'+str(i[3].idx), i.reschain_l, - [i.distance_aw, i.distance_dw], [i.d_angle, i.w_angle], i[0].coords, i[3].coords, - i.water.coords, i[7].residue.name+'_'+str(i[7].residue.idx)] - - Inter_list_all.append(Inter_list) - - for nr_k,k in enumerate(Inter_list_all): - LOGGER.info("%3i%22s%10s%26s%4s <---> %8s%12s%4s%6.1f" % (nr_k,k[0],k[1],k[2],k[3],k[4],k[5],k[6],k[7])) - - return Inter_list_all - - def showProteinInteractions_VMD(atoms, interactions, color='red',**kwargs): """Save information about protein interactions to a TCL file (filename) which can be further use in VMD to display all intercations in a graphical interface From d77f8b2b7a563b2466971a9e2ac0fee76bfb13c5 Mon Sep 17 00:00:00 2001 From: karolamik13 Date: Sun, 15 Oct 2023 21:31:34 +0200 Subject: [PATCH 04/11] fixed name of ligand in showLigandInteractions --- prody/proteins/interactions.py | 35 +++++++++++++++++----------------- 1 file changed, 18 insertions(+), 17 deletions(-) diff --git a/prody/proteins/interactions.py b/prody/proteins/interactions.py index 0c5f3f032..7dd7bf112 100644 --- a/prody/proteins/interactions.py +++ b/prody/proteins/interactions.py @@ -2029,41 +2029,41 @@ def showLigandInteractions(PLIP_output): param_inter = [method for method in dir(i) if method.startswith('_') is False] if str(type(i)).split('.')[-1].strip("'>") == 'hbond': - Inter_list = ['hbond',i.restype+str(i.resnr), i[0].type+'_'+str(i.d_orig_idx), i.reschain, - i.restype+str(i.resnr_l), i[2].type+'_'+str(i.a_orig_idx), i.reschain_l, + Inter_list = ['HBs',i.restype+str(i.resnr), i[0].type+'_'+str(i.d_orig_idx), i.reschain, + i.restype_l+str(i.resnr_l), i[2].type+'_'+str(i.a_orig_idx), i.reschain_l, i.distance_ad, i.angle, i[0].coords, i[2].coords] if str(type(i)).split('.')[-1].strip("'>") == 'saltbridge': - Inter_list = ['saltbridge',i.restype+str(i.resnr), '_'.join(map(str,i.negative.atoms_orig_idx)), i.reschain, - i.restype+str(i.resnr_l), '_'.join(map(str,i.positive.atoms_orig_idx)), i.reschain_l, + Inter_list = ['SBs',i.restype+str(i.resnr), '_'.join(map(str,i.negative.atoms_orig_idx)), i.reschain, + i.restype_l+str(i.resnr_l), '_'.join(map(str,i.positive.atoms_orig_idx)), i.reschain_l, i.distance, None, i.negative.center, i.positive.center] if str(type(i)).split('.')[-1].strip("'>") == 'pistack': - Inter_list = ['pistack',i.restype+str(i.resnr), '_'.join(map(str,i[0].atoms_orig_idx)), i.reschain, - i.restype+str(i.resnr_l), '_'.join(map(str,i[1].atoms_orig_idx)), i.reschain_l, + Inter_list = ['PiStack',i.restype+str(i.resnr), '_'.join(map(str,i[0].atoms_orig_idx)), i.reschain, + i.restype_l+str(i.resnr_l), '_'.join(map(str,i[1].atoms_orig_idx)), i.reschain_l, i.distance, i.angle, i[0].center, i[1].center] if str(type(i)).split('.')[-1].strip("'>") == 'pication': - Inter_list = ['pication',i.restype+str(i.resnr), '_'.join(map(str,i[0].atoms_orig_idx)), i.reschain, - i.restype+str(i.resnr_l), '_'.join(map(str,i[1].atoms_orig_idx)), i.reschain_l, + Inter_list = ['PiCat',i.restype+str(i.resnr), '_'.join(map(str,i[0].atoms_orig_idx)), i.reschain, + i.restype_l+str(i.resnr_l), '_'.join(map(str,i[1].atoms_orig_idx)), i.reschain_l, i.distance, None, i[0].center, i[1].center] if str(type(i)).split('.')[-1].strip("'>") == 'hydroph_interaction': - Inter_list = ['hydroph_interaction',i.restype+str(i.resnr), i[0].type+'_'+str(i[0].idx), i.reschain, - i.restype+str(i.resnr_l), i[2].type+'_'+str(i[2].idx), i.reschain_l, + Inter_list = ['HPh',i.restype+str(i.resnr), i[0].type+'_'+str(i[0].idx), i.reschain, + i.restype_l+str(i.resnr_l), i[2].type+'_'+str(i[2].idx), i.reschain_l, i.distance, None, i[0].coords, i[2].coords] if str(type(i)).split('.')[-1].strip("'>") == 'waterbridge': water = i.water - Inter_list = ['waterbridge',i.restype+str(i.resnr), i[0].type+'_'+str(i[0].idx), i.reschain, - i.restype+str(i.resnr_l), i[3].type+'_'+str(i[3].idx), i.reschain_l, + Inter_list = ['watBridge',i.restype+str(i.resnr), i[0].type+'_'+str(i[0].idx), i.reschain, + i.restype_l+str(i.resnr_l), i[3].type+'_'+str(i[3].idx), i.reschain_l, [i.distance_aw, i.distance_dw], [i.d_angle, i.w_angle], i[0].coords, i[3].coords, i.water.coords, i[7].residue.name+'_'+str(i[7].residue.idx)] Inter_list_all.append(Inter_list) for nr_k,k in enumerate(Inter_list_all): - LOGGER.info("%3i%22s%10s%26s%4s <---> %8s%12s%4s%6.1f" % (nr_k+1,k[0],k[1],k[2],k[3],k[4],k[5],k[6],k[7])) + LOGGER.info("%3i%12s%10s%26s%4s <---> %8s%12s%4s%6.1f" % (nr_k+1,k[0],k[1],k[2],k[3],k[4],k[5],k[6],k[7])) return Inter_list_all @@ -2130,15 +2130,15 @@ def calcLigandInteractions(atoms, **kwargs): if 'sele' in kwargs: sele = kwargs['sele'] else: - select='all not (water or protein or ion)' + sele='all not (water or protein or ion)' if 'ignore_ligs' in kwargs: ignore_ligs = kwargs['ignore_ligs'] else: ignore_ligs=['MAN'] - select = select+' and not (resname '+' '.join(ignore_ligs)+')' - ligand_select = atoms.select(select) + sele = sele+' and not (resname '+' '.join(ignore_ligs)+')' + ligand_select = atoms.select(sele) analyzedLigand = [] LOGGER.info("Detected ligands: ") for i in range(len(ligand_select.getResnums())): # It has to be done by each atom @@ -2155,7 +2155,8 @@ def calcLigandInteractions(atoms, **kwargs): Ligands.append(my_interactions) showLigandInteractions(my_interactions) except: - LOGGER.info(my_bsid+" not analyzed") + pass + #LOGGER.info(my_bsid+" not analyzed") return Ligands, analyzedLigand From e03b9c1d5f7e64ef5485486d82a3bc56b2f2dcbb Mon Sep 17 00:00:00 2001 From: karolamik13 Date: Mon, 16 Oct 2023 10:08:59 +0200 Subject: [PATCH 05/11] fixed bug with listing water bridges in ligand interactions --- prody/proteins/interactions.py | 35 ++++++++++++++++++++-------------- 1 file changed, 21 insertions(+), 14 deletions(-) diff --git a/prody/proteins/interactions.py b/prody/proteins/interactions.py index 7dd7bf112..a29801dad 100644 --- a/prody/proteins/interactions.py +++ b/prody/proteins/interactions.py @@ -2063,7 +2063,11 @@ def showLigandInteractions(PLIP_output): Inter_list_all.append(Inter_list) for nr_k,k in enumerate(Inter_list_all): - LOGGER.info("%3i%12s%10s%26s%4s <---> %8s%12s%4s%6.1f" % (nr_k+1,k[0],k[1],k[2],k[3],k[4],k[5],k[6],k[7])) + if k[0] == 'watBridge': + LOGGER.info("%3i%12s%10s%26s%4s <---> %8s%12s%4s%12s%14s" % (nr_k+1,k[0],k[1],k[2],k[3],k[4],k[5],k[6], + ' '.join(str(np.round(x, 2)) for x in k[7]), ' '.join(str(np.round(x, 2)) for x in k[8]))) + else: + LOGGER.info("%3i%12s%10s%26s%4s <---> %8s%12s%4s%6.1f" % (nr_k+1,k[0],k[1],k[2],k[3],k[4],k[5],k[6],k[7])) return Inter_list_all @@ -2135,30 +2139,30 @@ def calcLigandInteractions(atoms, **kwargs): if 'ignore_ligs' in kwargs: ignore_ligs = kwargs['ignore_ligs'] else: - ignore_ligs=['MAN'] + ignore_ligs=['MAN', 'SOD', 'CLA'] sele = sele+' and not (resname '+' '.join(ignore_ligs)+')' ligand_select = atoms.select(sele) analyzedLigand = [] - LOGGER.info("Detected ligands: ") - for i in range(len(ligand_select.getResnums())): # It has to be done by each atom - try: + + try: + for i in range(len(ligand_select.getResnums())): ResID = ligand_select.getResnames()[i] ChainID = ligand_select.getChids()[i] ResNames = ligand_select.getResnums()[i] my_bsid = str(ResID)+':'+str(ChainID)+':'+str(ResNames) - if my_bsid not in analyzedLigand: - LOGGER.info(my_bsid) + if my_bsid not in analyzedLigand: + LOGGER.info("LIGAND: {0}".format(my_bsid)) analyzedLigand.append(my_bsid) my_mol.analyze() my_interactions = my_mol.interaction_sets[my_bsid] # Contains all interaction data Ligands.append(my_interactions) showLigandInteractions(my_interactions) - except: - pass - #LOGGER.info(my_bsid+" not analyzed") - - return Ligands, analyzedLigand + + return Ligands, analyzedLigand + + except: + LOGGER.info("Ligand not found.") except ImportError: raise ImportError("Install Openbabel and PLIP.") @@ -2306,11 +2310,14 @@ def showLigandInteraction_VMD(atoms, interactions, **kwargs): except: filename = atoms.getTitle()+'_interaction.tcl' + pdb_name = atoms.getTitle()+'_sele.pdb' + writePDB(pdb_name, atoms) + tcl_file = open(filename, 'w') if len(interactions[0]) >= 10: - dic_color = {'hbond':'blue','pistack':'green','saltbridge':'yellow','pication':'orange', - 'hydroph_interaction':'silver','waterbridge':'cyan'} + dic_color = {'HBs':'blue','PiStack':'green','SBs':'yellow','PiCat':'orange', + 'HPh':'silver','watBridge':'cyan'} for i in interactions: tcl_file.write('draw color '+dic_color[i[0]]+'\n') From 34f4fb8e52495792e0292c70adef6db2c2bfdd72 Mon Sep 17 00:00:00 2001 From: karolamik13 Date: Mon, 16 Oct 2023 11:11:40 +0200 Subject: [PATCH 06/11] additional changes in protein-ligand interaction modules --- prody/proteins/interactions.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/prody/proteins/interactions.py b/prody/proteins/interactions.py index a29801dad..d0247dadd 100644 --- a/prody/proteins/interactions.py +++ b/prody/proteins/interactions.py @@ -2062,6 +2062,7 @@ def showLigandInteractions(PLIP_output): Inter_list_all.append(Inter_list) + LOGGER.info("%3s%12s%10s%20s%8s <---> %6s%10s%6s%10s%16s" % ('#','Type','Residue','Atoms','Chain','','Ligand','Atoms','Chain','Distance/Angle')) for nr_k,k in enumerate(Inter_list_all): if k[0] == 'watBridge': LOGGER.info("%3i%12s%10s%26s%4s <---> %8s%12s%4s%12s%14s" % (nr_k+1,k[0],k[1],k[2],k[3],k[4],k[5],k[6], @@ -2161,7 +2162,7 @@ def calcLigandInteractions(atoms, **kwargs): return Ligands, analyzedLigand - except: + except: LOGGER.info("Ligand not found.") except ImportError: From 0021ca3fce0768d5262ee7c2c476a99cdde0b9cf Mon Sep 17 00:00:00 2001 From: karolamik13 Date: Mon, 16 Oct 2023 12:27:16 +0200 Subject: [PATCH 07/11] LigandInteractionsTrajectory class --- prody/proteins/interactions.py | 106 ++++++++++++++++++++++++++++++++- 1 file changed, 105 insertions(+), 1 deletion(-) diff --git a/prody/proteins/interactions.py b/prody/proteins/interactions.py index d0247dadd..08d59e5c6 100644 --- a/prody/proteins/interactions.py +++ b/prody/proteins/interactions.py @@ -44,7 +44,7 @@ 'calcLigandInteractions', 'showLigandInteractions', 'showProteinInteractions_VMD', 'showLigandInteraction_VMD', 'calcHydrogenBondsTrajectory', 'calcHydrophobicOverlapingAreas', - 'Interactions', 'InteractionsTrajectory'] + 'Interactions', 'InteractionsTrajectory', 'LigandInteractionsTrajectory'] def cleanNumbers(listContacts): @@ -3627,4 +3627,108 @@ def getTimeInteractions(self, filename=None, **kwargs): return HBs, SBs, HPh, PiStack, PiCat, HPh, DiBs + + +class LigandInteractionsTrajectory(object): + + """Class for protein-ligand interaction analysis of DCD trajectory or multi-model PDB (Ensemble PDB).""" + + def __init__(self, name='Unknown'): + + self._atoms = None + self._traj = None + self._interactions_traj = None + + def calcLigandInteractionsTrajectory(self, atoms, trajectory=None, filename=None, **kwargs): + """Compute protein-ligand interactions for DCD trajectory or multi-model PDB + using PLIP library. + + :arg atoms: an Atomic object from which residues are selected + :type atoms: :class:`.Atomic` + + :arg trajectory: trajectory file + :type trajectory: class:`.Trajectory` + + :arg filename: Name of pkl filename in which interactions will be storage + :type filename: pkl + + :arg start_frame: index of first frame to read + :type start_frame: int + + :arg stop_frame: index of last frame to read + :type stop_frame: int """ + + + try: + coords = (atoms._getCoords() if hasattr(atoms, '_getCoords') else + atoms.getCoords()) + except AttributeError: + try: + checkCoords(coords) + except TypeError: + raise TypeError('coords must be an object ' + 'with `getCoords` method') + + interactions_all = [] + interactions_all_nb = [] + + start_frame = kwargs.pop('start_frame', 0) + stop_frame = kwargs.pop('stop_frame', -1) + + if trajectory is not None: + if isinstance(trajectory, Atomic): + trajectory = Ensemble(trajectory) + + nfi = trajectory._nfi + trajectory.reset() + numFrames = trajectory._n_csets + + if stop_frame == -1: + traj = trajectory[start_frame:] + else: + traj = trajectory[start_frame:stop_frame+1] + + atoms_copy = atoms.copy() + + for j0, frame0 in enumerate(traj, start=start_frame): + LOGGER.info(' ') + LOGGER.info('Frame: {0}'.format(j0)) + atoms_copy.setCoords(frame0.getCoords()) + + ligand_interactions, ligand = calcLigandInteractions(atoms_copy) + for ligs in range(len(ligand_interactions)-1): + interactions = showLigandInteractions(ligand_interactions[ligs]) + interactions_all.append(interactions) + interactions_all_nb.append(len(interactions)) + + else: + if atoms.numCoordsets() > 1: + for i in range(len(atoms.getCoordsets()[start_frame:stop_frame])): + LOGGER.info(' ') + LOGGER.info('Model: {0}'.format(i+start_frame)) + atoms.setACSIndex(i+start_frame) + + ligand_interactions, ligand = calcLigandInteractions(atoms) + for ligs in range(len(ligand_interactions)-1): + interactions = showLigandInteractions(ligand_interactions[ligs]) + interactions_all.append(interactions) + interactions_all_nb.append(len(interactions)) + + else: + LOGGER.info('Include trajectory or use multi-model PDB file.') + + self._atoms = atoms + self._traj = trajectory + self._interactions_traj = interactions_all + self._interactions_nb_traj = interactions_all_nb + + if filename is not None: + import pickle + with open(str(filename)+'.pkl', 'wb') as f: + pickle.dump(self._interactions_traj, f) + LOGGER.info('File with interactions saved.') + + return interactions_all + + From 9489bffa0ea817da0abeb8020abd9e95e1eb1b13 Mon Sep 17 00:00:00 2001 From: karolamik13 Date: Mon, 16 Oct 2023 21:46:22 +0200 Subject: [PATCH 08/11] New class LigandInteractionsTrajectory is created --- prody/proteins/interactions.py | 45 ++++++++++++++++++++++++++-------- 1 file changed, 35 insertions(+), 10 deletions(-) diff --git a/prody/proteins/interactions.py b/prody/proteins/interactions.py index 08d59e5c6..b01176f0b 100644 --- a/prody/proteins/interactions.py +++ b/prody/proteins/interactions.py @@ -3690,18 +3690,22 @@ def calcLigandInteractionsTrajectory(self, atoms, trajectory=None, filename=None traj = trajectory[start_frame:stop_frame+1] atoms_copy = atoms.copy() - + for j0, frame0 in enumerate(traj, start=start_frame): LOGGER.info(' ') LOGGER.info('Frame: {0}'.format(j0)) atoms_copy.setCoords(frame0.getCoords()) ligand_interactions, ligand = calcLigandInteractions(atoms_copy) - for ligs in range(len(ligand_interactions)-1): - interactions = showLigandInteractions(ligand_interactions[ligs]) - interactions_all.append(interactions) - interactions_all_nb.append(len(interactions)) - + + ligs_per_frame_interactions = [] + for ligs in ligand_interactions: + LP_interactions = showLigandInteractions(ligs) + ligs_per_frame_interactions.extend(LP_interactions) + + interactions_all.append(ligs_per_frame_interactions) + interactions_all_nb.append(len(ligs_per_frame_interactions)) + else: if atoms.numCoordsets() > 1: for i in range(len(atoms.getCoordsets()[start_frame:stop_frame])): @@ -3710,11 +3714,15 @@ def calcLigandInteractionsTrajectory(self, atoms, trajectory=None, filename=None atoms.setACSIndex(i+start_frame) ligand_interactions, ligand = calcLigandInteractions(atoms) - for ligs in range(len(ligand_interactions)-1): - interactions = showLigandInteractions(ligand_interactions[ligs]) - interactions_all.append(interactions) - interactions_all_nb.append(len(interactions)) + ligs_per_frame_interactions = [] + for ligs in ligand_interactions: + LP_interactions = showLigandInteractions(ligs) + ligs_per_frame_interactions.extend(LP_interactions) + + interactions_all.append(ligs_per_frame_interactions) + interactions_all_nb.append(len(ligs_per_frame_interactions)) + else: LOGGER.info('Include trajectory or use multi-model PDB file.') @@ -3732,3 +3740,20 @@ def calcLigandInteractionsTrajectory(self, atoms, trajectory=None, filename=None return interactions_all + def getLigandInteractions(self, **kwargs): + """Return the list of protein-ligand interactions.""" + + return self._interactions_traj + + + def getAtoms(self): + """Returns associated atoms.""" + + return self._atoms + + + def getInteractionsNumber(self): + """Return the number of interactions in each frame.""" + + return self._interactions_nb_traj + From ada74374f680b57101b6220cb75c5f12301af404 Mon Sep 17 00:00:00 2001 From: karolamik13 Date: Tue, 17 Oct 2023 08:51:19 +0200 Subject: [PATCH 09/11] parseLigandInteractions added to LigandInteractionsTrajectory class --- prody/proteins/interactions.py | 19 +++++++++++++++++-- 1 file changed, 17 insertions(+), 2 deletions(-) diff --git a/prody/proteins/interactions.py b/prody/proteins/interactions.py index b01176f0b..eba9e410a 100644 --- a/prody/proteins/interactions.py +++ b/prody/proteins/interactions.py @@ -3628,7 +3628,6 @@ def getTimeInteractions(self, filename=None, **kwargs): return HBs, SBs, HPh, PiStack, PiCat, HPh, DiBs - class LigandInteractionsTrajectory(object): @@ -3752,8 +3751,24 @@ def getAtoms(self): return self._atoms - def getInteractionsNumber(self): + def getLigandInteractionsNumber(self): """Return the number of interactions in each frame.""" return self._interactions_nb_traj + + def parseLigandInteractions(self, filename): + """Import interactions from analysis of trajectory which was saved via + calcLigandInteractionsTrajectory(). + + :arg filename: Name of pkl file in which interactions will be storage + :type filename: pkl """ + + import pickle + with open(filename, 'rb') as f: + data = pickle.load(f) + + self._interactions_traj = data + self._interactions_nb_traj = [[len(sublist) if sublist else 0 for sublist in sublist] for sublist in data] + + return data From aecfea2e0cedea6c08749d599b309ad6496a2660 Mon Sep 17 00:00:00 2001 From: karolamik13 Date: Tue, 17 Oct 2023 08:59:08 +0200 Subject: [PATCH 10/11] More info in LigandInteractionsTrajectory docx --- prody/proteins/interactions.py | 17 +++++++++++++++-- 1 file changed, 15 insertions(+), 2 deletions(-) diff --git a/prody/proteins/interactions.py b/prody/proteins/interactions.py index eba9e410a..e8969463d 100644 --- a/prody/proteins/interactions.py +++ b/prody/proteins/interactions.py @@ -3631,7 +3631,20 @@ def getTimeInteractions(self, filename=None, **kwargs): class LigandInteractionsTrajectory(object): - """Class for protein-ligand interaction analysis of DCD trajectory or multi-model PDB (Ensemble PDB).""" + """Class for protein-ligand interaction analysis of DCD trajectory or multi-model PDB (Ensemble PDB). + This class is using PLIP to provide the interactions. Install PLIP before using it. + + ## Instalation of PLIP using conda: + >>> conda install -c conda-forge plip + ## https://anaconda.org/conda-forge/plip + # https://github.com/pharmai/plip/blob/master/DOCUMENTATION.md + + ## Instalation using PIP: + >>> pip install plip + + .. [SS15] Salentin S., Schreiber S., Haupt V. J., Adasme M. F., Schroeder M. + PLIP: fully automated protein–ligand interaction profiler + *Nucl. Acids Res.* **2015** 43:W443-W447. """ def __init__(self, name='Unknown'): @@ -3657,7 +3670,6 @@ def calcLigandInteractionsTrajectory(self, atoms, trajectory=None, filename=None :arg stop_frame: index of last frame to read :type stop_frame: int """ - try: coords = (atoms._getCoords() if hasattr(atoms, '_getCoords') else @@ -3772,3 +3784,4 @@ def parseLigandInteractions(self, filename): self._interactions_nb_traj = [[len(sublist) if sublist else 0 for sublist in sublist] for sublist in data] return data + From ade09fd222b5cb7d41e10a3e01d875506d212d81 Mon Sep 17 00:00:00 2001 From: karolamik13 Date: Tue, 17 Oct 2023 20:33:33 +0200 Subject: [PATCH 11/11] changes in import plip and listLigandInteractions name --- prody/proteins/interactions.py | 103 +++++++++++++++++---------------- 1 file changed, 52 insertions(+), 51 deletions(-) diff --git a/prody/proteins/interactions.py b/prody/proteins/interactions.py index e8969463d..f3712d94e 100644 --- a/prody/proteins/interactions.py +++ b/prody/proteins/interactions.py @@ -41,7 +41,7 @@ 'calcPiCationTrajectory', 'calcHydrophobicTrajectory', 'calcDisulfideBondsTrajectory', 'calcProteinInteractions', 'calcStatisticsInteractions', 'calcDistribution', 'calcSASA', 'calcVolume','compareInteractions', 'showInteractionsGraph', - 'calcLigandInteractions', 'showLigandInteractions', + 'calcLigandInteractions', 'listLigandInteractions', 'showProteinInteractions_VMD', 'showLigandInteraction_VMD', 'calcHydrogenBondsTrajectory', 'calcHydrophobicOverlapingAreas', 'Interactions', 'InteractionsTrajectory', 'LigandInteractionsTrajectory'] @@ -2014,7 +2014,7 @@ def calcDistribution(interactions, residue1, residue2=None, **kwargs): LOGGER.info(i) -def showLigandInteractions(PLIP_output): +def listLigandInteractions(PLIP_output): """Create a list of interactions from PLIP output created using calcLigandInteractions(). Results can be displayed in VMD. @@ -2115,58 +2115,59 @@ def calcLigandInteractions(atoms, **kwargs): raise TypeError('coords must be an object ' 'with `getCoords` method') try: - from plip.structure.preparation import PDBComplex - import tempfile + from plip.structure.preparation import PDBComplex - with tempfile.NamedTemporaryFile(delete=False, suffix=".pdb") as temp_pdb_file: - writePDB(temp_pdb_file.name, atoms) - temp_pdb_file_name = temp_pdb_file.name + except ImportError: + raise ImportError("Install Openbabel and PLIP.") + + import tempfile + with tempfile.NamedTemporaryFile(delete=False, suffix=".pdb") as temp_pdb_file: + writePDB(temp_pdb_file.name, atoms) + temp_pdb_file_name = temp_pdb_file.name - try: - if atoms.hydrogen == None or atoms.hydrogen.numAtoms() < 30: # if there is no hydrogens in PDB structure - LOGGER.info("Lack of hydrogens in the structure. Use addMissingAtoms to add them.") - except: - pass + try: + if atoms.hydrogen == None or atoms.hydrogen.numAtoms() < 30: # if there is no hydrogens in PDB structure + LOGGER.info("Lack of hydrogens in the structure. Use addMissingAtoms to add them.") + except: + pass - Ligands = [] # Ligands can be more than one - my_mol = PDBComplex() - my_mol.load_pdb(temp_pdb_file_name) # Load the PDB file into PLIP class - - if 'sele' in kwargs: - sele = kwargs['sele'] - else: - sele='all not (water or protein or ion)' + Ligands = [] # Ligands can be more than one + my_mol = PDBComplex() + my_mol.load_pdb(temp_pdb_file_name) # Load the PDB file into PLIP class + + if 'sele' in kwargs: + sele = kwargs['sele'] + else: + sele='all not (water or protein or ion)' - if 'ignore_ligs' in kwargs: - ignore_ligs = kwargs['ignore_ligs'] - else: - ignore_ligs=['MAN', 'SOD', 'CLA'] - - sele = sele+' and not (resname '+' '.join(ignore_ligs)+')' - ligand_select = atoms.select(sele) - analyzedLigand = [] - - try: - for i in range(len(ligand_select.getResnums())): - ResID = ligand_select.getResnames()[i] - ChainID = ligand_select.getChids()[i] - ResNames = ligand_select.getResnums()[i] - my_bsid = str(ResID)+':'+str(ChainID)+':'+str(ResNames) - if my_bsid not in analyzedLigand: - LOGGER.info("LIGAND: {0}".format(my_bsid)) - analyzedLigand.append(my_bsid) - my_mol.analyze() - my_interactions = my_mol.interaction_sets[my_bsid] # Contains all interaction data - Ligands.append(my_interactions) - showLigandInteractions(my_interactions) - - return Ligands, analyzedLigand - - except: - LOGGER.info("Ligand not found.") + if 'ignore_ligs' in kwargs: + ignore_ligs = kwargs['ignore_ligs'] + else: + ignore_ligs=['MAN', 'SOD', 'CLA'] + + sele = sele+' and not (resname '+' '.join(ignore_ligs)+')' + ligand_select = atoms.select(sele) + analyzedLigand = [] + + try: + for i in range(len(ligand_select.getResnums())): + ResID = ligand_select.getResnames()[i] + ChainID = ligand_select.getChids()[i] + ResNames = ligand_select.getResnums()[i] + my_bsid = str(ResID)+':'+str(ChainID)+':'+str(ResNames) + if my_bsid not in analyzedLigand: + LOGGER.info("LIGAND: {0}".format(my_bsid)) + analyzedLigand.append(my_bsid) + my_mol.analyze() + my_interactions = my_mol.interaction_sets[my_bsid] # Contains all interaction data + Ligands.append(my_interactions) + listLigandInteractions(my_interactions) + + return Ligands, analyzedLigand + + except: + LOGGER.info("Ligand not found.") - except ImportError: - raise ImportError("Install Openbabel and PLIP.") def showProteinInteractions_VMD(atoms, interactions, color='red',**kwargs): @@ -3711,7 +3712,7 @@ def calcLigandInteractionsTrajectory(self, atoms, trajectory=None, filename=None ligs_per_frame_interactions = [] for ligs in ligand_interactions: - LP_interactions = showLigandInteractions(ligs) + LP_interactions = listLigandInteractions(ligs) ligs_per_frame_interactions.extend(LP_interactions) interactions_all.append(ligs_per_frame_interactions) @@ -3728,7 +3729,7 @@ def calcLigandInteractionsTrajectory(self, atoms, trajectory=None, filename=None ligs_per_frame_interactions = [] for ligs in ligand_interactions: - LP_interactions = showLigandInteractions(ligs) + LP_interactions = listLigandInteractions(ligs) ligs_per_frame_interactions.extend(LP_interactions) interactions_all.append(ligs_per_frame_interactions)