chemosim-lab / ProLIF

Interaction Fingerprints for protein-ligand complexes and more
https://prolif.readthedocs.io
Apache License 2.0
365 stars 68 forks source link

KeyError: 'PiStacking' #76

Closed amir-tagh closed 1 year ago

amir-tagh commented 2 years ago

Hello,

I am using the following steps to get the interactions between an RNA and a ligand and in getting the pi-stacking interactions I am getting an error message. Can you please let me know how to solve this.

Thanks, Amir

raise KeyError(key) KeyError: 'PiStacking'

!/usr/bin/env python

from rdkit import Chem import sys from rdkit import Chem from rdkit.Chem import Draw import MDAnalysis as mda import prolif as plf import pandas as pd import numpy as np from prolif.plotting.network import LigNetwork

load protein

prot = mda.Universe(plf.datafiles.datapath / "vina" / "rec.pdb")

prot = mda.Universe("RNA.pdb") prot = plf.Molecule.from_mda(prot) print(prot.n_residues)

load ligands

path = str("ligand.sdf") lig_suppl = plf.sdf_supplier(path)

generate fingerprint

fp = plf.Fingerprint() fp.run_from_iterable(lig_suppl, prot) df = fp.to_dataframe() print(df)

mol = Chem.MolFromPDBFile("RNA.pdb", removeHs=False) prot = plf.Molecule(mol) mol = Chem.MolFromPDBFile("ligand.pdb", removeHs=False) lig = plf.Molecule(mol) fp = plf.Fingerprint() fp.run_from_iterable([lig], prot,residues="all")

df = fp.to_dataframe()

show only the 10 first frames

print(df.head(10))

drop the ligand residue column since there's only a single ligand residue

df = df.droplevel("ligand", axis=1) print(df.head(5))

show all pi-stacking interactions

df.xs("PiStacking", level="interaction", axis=1).head(5)

cbouy commented 2 years ago

Hi Amir,

I'll need a bit more info to solve this. The error means that there is no PiStacking interaction detected, but it could come from multiple things:

from rdkit import Chem import sys from rdkit import Chem from rdkit.Chem import Draw import MDAnalysis as mda import prolif as plf import pandas as pd import numpy as np from prolif.plotting.network import LigNetwork

load protein

prot = mda.Universe("RNA.pdb") prot = plf.Molecule.from_mda(prot) print(prot.n_residues)

load ligands

path = str("ligand.sdf") lig_suppl = plf.sdf_supplier(path)

generate fingerprint

fp = plf.Fingerprint() fp.run_from_iterable(lig_suppl, prot) df = fp.to_dataframe() print(df.head()) print(df.xs("PiStacking", level="interaction", axis=1).head(5))


Also the second block WILL NOT WORK:
```python
mol = Chem.MolFromPDBFile("RNA.pdb", removeHs=False)
prot = plf.Molecule(mol)
mol = Chem.MolFromPDBFile("ligand.pdb", removeHs=False)
lig = plf.Molecule(mol)
fp = plf.Fingerprint()
fp.run_from_iterable([lig], prot,residues="all")
df = fp.to_dataframe()

Because RDKit cannot determine bond order information and formal charge (which ProLIF needs) from any ligand type, it only does this for the standard protein residues. To make it work, you should replace this:

mol = Chem.MolFromPDBFile("ligand.pdb", removeHs=False)
lig = plf.Molecule(mol)

with this:

lig = mda.Universe("ligand.pdb")
lig = plf.Molecule.from_mda(lig)

Best, Cedric

amir-tagh commented 2 years ago

Hello Cedric,

Thanks a lot for your help. I just trimmed done the code just to the second block. I should mention that I am running this on an "RNA" and not a protein, I am not sure if this applicable to RNA at all. I can see the stacking interactions with another tool, just to make sure stacking interactions are there but not with Prolif. the whole error message is copied.

Best, Amir

does your RNA.pdb file contain all hydrogen atoms explicitly? Yes same question for ligand.pdb Yes

!/usr/bin/env python

from rdkit import Chem import sys from rdkit import Chem from rdkit.Chem import Draw import MDAnalysis as mda import prolif as plf import pandas as pd import numpy as np from prolif.plotting.network import LigNetwork

mol = Chem.MolFromPDBFile("RNA.pdb", removeHs=False) prot = plf.Molecule(mol)

lig = mda.Universe("ligand.pdb") lig = plf.Molecule.from_mda(lig)

fp = plf.Fingerprint() fp.run_from_iterable([lig], prot,residues="all") df = fp.to_dataframe() print(df)

print(df.xs("PiStacking", level="interaction", axis=1).head(5))

'''ligand LIG25
protein RG5. RA6. RC7. RG19. RC20. interaction Hydrophobic Hydrophobic HBDonor Hydrophobic Hydrophobic Hydrophobic Frame
0 True True True True True True Traceback (most recent call last): File "/home/amir/Projects/P4D6_docking/fingeRNA_dock10/Prolif_analysis/./Prolif.py", line 43, in print(df.xs("PiStacking", level="interaction", axis=1).head(5)) File "/home/amir/anaconda3/envs/prolif/lib/python3.10/site-packages/pandas/core/generic.py", line 3836, in xs loc, new_ax = labels.get_loc_level(key, level=level, drop_level=drop_level) File "/home/amir/anaconda3/envs/prolif/lib/python3.10/site-packages/pandas/core/indexes/multi.py", line 2972, in get_loc_level loc, mi = self._get_loc_level(key, level=level) File "/home/amir/anaconda3/envs/prolif/lib/python3.10/site-packages/pandas/core/indexes/multi.py", line 3113, in _get_loc_level indexer = self._get_level_indexer(key, level=level) File "/home/amir/anaconda3/envs/prolif/lib/python3.10/site-packages/pandas/core/indexes/multi.py", line 3235, in _get_level_indexer raise KeyError(key) KeyError: 'PiStacking' '''

cbouy commented 2 years ago

Since both files have all hydrogen atoms explicit and there shouldn't be anything wrong with the code snippet, I think it might be the default parameters for that interaction which aren't suitable in your case. You could try tweaking a bit the parameters for that specific interaction by adding the following snippet after importing the libraries:

class PiStacking(plf.interactions.PiStacking):
    def __init__(self):
        super().__init__(centroid_distance=6.2, shortest_distance=4.0)

This will automatically overwrite the parameters for the PiStacking interaction. If that's still not enough, try slightly increasing the 2 distance parameters by 0.2 again.

And if this still doesn't work, it might the current implementation of Pi-Stacking detection which is misbehaving in that particular case. I have plans to change it at some point but there's nothing you can do for now unfortunately.

Best, Cedric

amir-tagh commented 2 years ago

Thanks, I used the snippet for overwriting the parameters but couldn't get the stacking interactions. I will try changing the distance parameters see if I can get it. Thanks for your help.

best, Amir