diff --git a/modules/mol/alg/pymod/ligand_scoring.py b/modules/mol/alg/pymod/ligand_scoring.py index 6a9895d6b6addf5c2456c726259a12c7b3f05ef7..d1d1f2622decf98cbc72d4e2239b1f62311d288c 100644 --- a/modules/mol/alg/pymod/ligand_scoring.py +++ b/modules/mol/alg/pymod/ligand_scoring.py @@ -317,11 +317,13 @@ class LigandScorer: self._unassigned_model_ligand_short = None self._unassigned_target_ligand_descriptions = None self._unassigned_model_ligand_descriptions = None - # Keep track of symmetries/isomorphisms + # Keep track of symmetries/isomorphisms (regardless of scoring) # 0.0: no isomorphism # 1.0: isomorphic - # np.nan: not assessed yet + # np.nan: not assessed yet - that's why we can't use a boolean self._assignment_isomorphisms = None + # Keep track of match coverage (only in case there was a score) + self._assignment_match_coverage = None if custom_mapping is not None: self._set_custom_mapping(custom_mapping) @@ -620,6 +622,8 @@ class LigandScorer: (len(self.target_ligands), len(self.model_ligands)), dtype=dict) self._assignment_isomorphisms = np.full( (len(self.target_ligands), len(self.model_ligands)), fill_value=np.nan) + self._assignment_match_coverage = np.full( + (len(self.target_ligands), len(self.model_ligands)), fill_value=np.nan) for target_i, target_ligand in enumerate(self.target_ligands): LogVerbose("Analyzing target ligand %s" % target_ligand) @@ -649,7 +653,6 @@ class LigandScorer: by_atom_index=True) LogVerbose("Ligands %s and %s symmetry match" % ( str(model_ligand), str(target_ligand))) - self._assignment_isomorphisms[target_i, model_i] = 1 except NoSymmetryError: # Ligands are different - skip LogVerbose("No symmetry between %s and %s" % ( @@ -661,6 +664,9 @@ class LigandScorer: continue substructure_match = len(symmetries[0][0]) != len( model_ligand.atoms) + coverage = len(symmetries[0][0]) / len(model_ligand.atoms) + self._assignment_match_coverage = coverage + self._assignment_isomorphisms[target_i, model_i] = 1 rmsd = SCRMSD(model_ligand, target_ligand, transformation=binding_site.transform, @@ -682,6 +688,7 @@ class LigandScorer: "chain_mapping": binding_site.GetFlatChainMapping(), "transform": binding_site.transform, "substructure_match": substructure_match, + "coverage": coverage, "inconsistent_residues": binding_site.inconsistent_residues, } if self.unassigned: @@ -746,6 +753,7 @@ class LigandScorer: "chain_mapping": binding_site.GetFlatChainMapping(), "transform": binding_site.transform, "substructure_match": substructure_match, + "coverage": coverage, "inconsistent_residues": binding_site.inconsistent_residues, } if self.unassigned: @@ -1130,6 +1138,9 @@ class LigandScorer: (substructure) match. A value of `True` indicates that the target ligand covers only part of the model, while `False` indicates a perfect match. + * `coverage`: the fraction of model atoms covered by the assigned + target ligand, in the interval (0, 1]. If `substructure_match` + is `False`, this will always be 1. * `inconsistent_residues`: a list of tuples of mapped residues views (:class:`~ost.mol.ResidueView`) with residue names that differ between the reference and the model, respectively. diff --git a/modules/mol/alg/tests/test_ligand_scoring.py b/modules/mol/alg/tests/test_ligand_scoring.py index 65249c2b54d47af75777d725ce9d7621df203088..51fb241dc3411ac9841d2693df52a32654bdd40b 100644 --- a/modules/mol/alg/tests/test_ligand_scoring.py +++ b/modules/mol/alg/tests/test_ligand_scoring.py @@ -596,6 +596,26 @@ class TestLigandScoring(unittest.TestCase): unassigned=True, rmsd_assignment=True) + def test_substructure_match(self): + """Test that substructure_match=True works.""" + trg = _LoadMMCIF("1r8q.cif.gz") + mdl = _LoadMMCIF("P84080_model_02.cif.gz") + + trg_g3d1 = trg.FindResidue("F", 1) + mdl_g3d = mdl.FindResidue("L_2", 1) + + # Skip PA, PB and O[1-3]A and O[1-3]B in target and model + # ie 8 / 32 atoms => coverage 0.75 + # We assume atom index are fixed and won't change + trg_g3d1_sub_ent = trg_g3d1.Select("aindex>6019") + trg_g3d1_sub = trg_g3d1_sub_ent.residues[0] + + # Substructure matches + sc = LigandScorer(mdl.Select("protein=True"), trg.Select("protein=True"), + model_ligands=[mdl_g3d], target_ligands=[trg_g3d1_sub], + substructure_match=True) + assert sc.rmsd_details["L_2"][1]["coverage"] == 0.75 + if __name__ == "__main__": from ost import testutils if testutils.DefaultCompoundLibIsSet():