chemosim-lab / ProLIF

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

HBA/HBD not showing #144

Closed JenkeScheen closed 1 year ago

JenkeScheen commented 1 year ago

Hi,

I'm having trouble getting ProLIF to detect some HBAs/HBDs in a structure. With PyMol using the show_contacts.py plugin: image

Quite a bit going on there. However, running with ProLIF:

prot = plf.Molecule(Chem.MolFromPDBFile("mac1_prot_hpp.pdb", removeHs=False))
lig = plf.Molecule(Chem.MolFromPDBFile("mac1_lig_truncated.pdb", removeHs=False))

prolif_fp = plf.Fingerprint()
prolif_fp.run_from_iterable([lig], prot, n_jobs=1)

# print out every interaction type detected by ProLIF:
[ print(list(data.keys())[0]) for data in list(prolif_fp.ifp[0].values()) ]

returns:

Hydrophobic
Hydrophobic
Hydrophobic
VdWContact
Hydrophobic

I tried playing around with plf.Fingerprint(parameters={"HBAcceptor": {"distance": x, "DHA_angle":(y, z)}}) but no settings (even extreme ones) enable detection of the interactions. Am I missing something here? Is there an issue with the protonation states that show_contacts.py is dealing with but ProLIF isn't?

Input files: issue_reproduce_inputs.tar.gz

JenkeScheen commented 1 year ago

this is with prolif 2.0.0.dev0 pypi_0 pypi

MIHassanLab commented 1 year ago

While using PDB files for interaction fingerprint and subsequent LigNetwork analysis, I am also facing the same issue. Only the VdWContacts and Hydrophobic interactions are returned.

cbouy commented 1 year ago

Hi,

In this particular case the protein is correctly prepared (i.e. all hydrogens explicit) but not the ligand so RDKit reads it in an unexpected way (incorrect hydrogen assignment circled): image instead of what I believe should be this image

If you have a SMILES string of the original compound you can fix the problem like so:

from rdkit.Chem.AllChem import AssignBondOrdersFromTemplate

reference = Chem.MolFromSmiles("CNC(=O)c1cc2cncnc2[nH]1")
rdlig = Chem.MolFromPDBFile("mac1_lig_truncated.pdb", removeHs=False)

mol = AssignBondOrdersFromTemplate(reference, rdlig)
mol = Chem.RemoveHs(mol)
Chem.SanitizeMol(mol)
mol = Chem.AddHs(mol, addCoords=True, addResidueInfo=True)
lig = plf.Molecule.from_rdkit(mol)

and this will detect the other interactions correctly.

In general I'd avoid using PDB format for the ligand if possible, it's really not made for small molecules...

Bonus: Since you have installed the development version for the next release (2.0.0) which includes the detection of all interactions occurrences (install of a single one per residue) you can do the following to display both interactions with that ASP residue:

from prolif.plotting.network import LigNetwork

fp = plf.Fingerprint(count=True)
fp.run_from_iterable([lig], prot, n_jobs=1)
net = LigNetwork.from_fingerprint(fp, lig, kind="frame", display_all=True)
net.display()

which should display this (only showing the HBond interactions for clarity) image

Without count=True you would only see one HBond for the ASP.

Hope this clarifies things, please close the issue if you don't have any other related questions.

Best, Cédric

JenkeScheen commented 1 year ago

Thanks a lot for this, I'll make sure to get protonation right for these.

In general I'd avoid using PDB format for the ligand if possible, it's really not made for small molecules... Agreed, I was getting SegFaults while using SDFs though, perhaps something is off with the input SDF.