ProLIF icon indicating copy to clipboard operation
ProLIF copied to clipboard

Interactions with a Mg atom

Open martutebon opened this issue 3 years ago • 1 comments

Hi, Thank you very much for your work! Very useful piece of Software.

I am experiencing the following issue. In my system there's a magnesium atom and the ligand is interacting with it. The selections are working.

lig = u.select_atoms("resname LIG") mg = u.select_atoms("resname MG")

Nevertheless, when I try to analyze the fingergerprints, I get: AttributeError: No hydrogen atom could be found in the topology, but the converter requires all hydrogens to be explicit. You can use the parameter NoImplicit=False when using the converter to allow implicit hydrogens and disable inferring bond orders and charges. You can also use force=True to ignore this error.

ligand/protein and ligand/water interactions are working fine.

Am I doing something wrong?

Thank you very much!

martutebon avatar Jun 11 '22 22:06 martutebon

Hi @martutebon and thank you!

Am I doing something wrong?

Nope! It's just that by default the RDKitConverter in MDAnalysis was set to automatically fail if a user is trying to convert a selection that doesn't have any hydrogen atom, since the RDKitConverter requires all hydrogen atoms to be explicit. Since in your case it's normal that the Magnesium selection doesn't have any H atom, you'd have to pass force=True to the RDKitConverter.

I'm assuming you're doing the fingerprint analysis on a trajectory. Unfortunately the functionality to pass arguments to the RDKitConverter is not directly exposed when using fp.run, it's only available when you manually create a molecule with plf.Molecule.from_mda (which fp.run uses internally). However here's a little "hack" to replace the plf.Molecule.from_mda function with one that automatically passes the force=True argument to the RDKitConverter:


import prolif as plf
import MDAnalysis as mda

def custom_from_mda(obj, selection=None, **kwargs):
     ag = obj.select_atoms(selection) if selection else obj.atoms
     assert ag.n_atoms > 0, "Empty atomgroup"
     mol = ag.convert_to.rdkit(force=True, **kwargs)
     return plf.Molecule(mol)

plf.Molecule.from_mda = custom_from_mda

# rests of your analysis script after this line

Then you should be able to use ProLIF as usual. I'll try to add the possibility to pass arguments to the converter directly from fp.run in the future.

Tell me if that works for you.

cbouy avatar Jun 13 '22 20:06 cbouy