6 from ost
import settings
8 from ost
import LogScript, LogInfo, LogDebug
20 from ost
import bindings
26 """Scorer specific for a reference/model pair
28 Finds best possible binding site representation of reference in model given
29 lDDT score. Uses :class:`ost.mol.alg.chain_mapping.ChainMapper` to deal with
32 :param reference: Reference structure
33 :type reference: :class:`ost.mol.EntityView`/:class:`ost.mol.EntityHandle`
34 :param model: Model structure
35 :type model: :class:`ost.mol.EntityView`/:class:`ost.mol.EntityHandle`
36 :param residue_number_alignment: Passed to ChainMapper constructor
37 :type residue_number_alignment: :class:`bool`
40 residue_number_alignment=False):
42 resnum_alignments=residue_number_alignment)
46 def ScoreBS(self, ligand, radius = 4.0, lddt_radius=10.0):
47 """Computes binding site lDDT score given *ligand*. Best possible
48 binding site representation is selected by lDDT but other scores such as
49 CA based RMSD and GDT are computed too and returned.
51 :param ligand: Defines the scored binding site, i.e. provides positions
52 to perform proximity search
53 :type ligand: r'((Residue)|(Chain)|(Entity))((View)|(Handle))'
54 :param radius: Reference residues with any atom position within *radius*
55 of *ligand* consitute the scored binding site
56 :type radius: :class:`float`
57 :param lddt_radius: Passed as *inclusion_radius* to
58 :class:`ost.mol.alg.lddt.lDDTScorer`
59 :type lddt_radius: :class:`float`
60 :returns: Object of type :class:`ost.mol.alg.chain_mapping.ReprResult`
61 containing all atom lDDT score and mapping information.
62 None if no representation could be found.
66 ref_residues_hashes = set()
67 for ligand_at
in ligand.atoms:
68 close_atoms = self.
refref.FindWithin(ligand_at.GetPos(), radius)
69 for close_at
in close_atoms:
70 ref_res = close_at.GetResidue()
71 h = ref_res.handle.GetHashCode()
72 if h
not in ref_residues_hashes:
73 ref_residues_hashes.add(h)
78 ref_bs = self.
refref.CreateEmptyView()
79 for ch
in self.
refref.chains:
81 if r.handle.GetHashCode()
in ref_residues_hashes:
82 ref_bs.AddResidue(r, mol.ViewAddFlag.INCLUDE_ALL)
86 inclusion_radius = lddt_radius)
94 """ Helper class to access the various scores available from ost.mol.alg
96 Deals with structure cleanup, chain mapping, interface identification etc.
97 Intermediate results are available as attributes.
99 :param model: Model structure - a deep copy is available as :attr:`~model`.
100 Additionally, :func:`ost.mol.alg.Molck` using *molck_settings*
102 :type model: :class:`ost.mol.EntityHandle`/:class:`ost.mol.EntityView`
103 :param target: Target structure - a deep copy is available as :attr:`~target`.
104 Additionally, :func:`ost.mol.alg.Molck` using *molck_settings*
106 :type target: :class:`ost.mol.EntityHandle`/:class:`ost.mol.EntityView`
107 :param resnum_alignments: Whether alignments between chemically equivalent
108 chains in *model* and *target* can be computed
109 based on residue numbers. This can be assumed in
110 benchmarking setups such as CAMEO/CASP.
111 :type resnum_alignments: :class:`bool`
112 :param molck_settings: Settings used for Molck on *model* and *target*, if
113 set to None, a default object is constructed by
114 setting everything except rm_zero_occ_atoms and
116 :class:`ost.mol.alg.MolckSettings` constructor.
117 :type molck_settings: :class:`ost.mol.alg.MolckSettings`
118 :param cad_score_exec: Explicit path to voronota-cadscore executable from
119 voronota installation from
120 https://github.com/kliment-olechnovic/voronota. If
121 not given, voronota-cadscore must be in PATH if any
122 of the CAD score related attributes is requested.
123 :type cad_score_exec: :class:`str`
124 :param custom_mapping: Provide custom chain mapping between *model* and
125 *target*. Dictionary with target chain names as key
126 and model chain names as value.
127 :attr:`~mapping` is constructed from this.
128 :type custom_mapping: :class:`dict`
129 :param custom_rigid_mapping: Provide custom chain mapping between *model*
130 and *target*. Dictionary with target chain
131 names as key and model chain names as value.
132 :attr:`~rigid_mapping` is constructed from
134 :type custom_rigid_mapping: :class:`dict`
135 :param usalign_exec: Explicit path to USalign executable used to compute
136 TM-score. If not given, TM-score will be computed
137 with OpenStructure internal copy of USalign code.
138 :type usalign_exec: :class:`str`
139 :param lddt_no_stereochecks: Whether to compute lDDT without stereochemistry
141 :type lddt_no_stereochecks: :class:`bool`
142 :param n_max_naive: Parameter for chain mapping. If the number of possible
143 mappings is <= *n_max_naive*, the full
144 mapping solution space is enumerated to find the
145 the optimum. A heuristic is used otherwise. The default
146 of 40320 corresponds to an octamer (8! = 40320).
147 A structure with stoichiometry A6B2 would be
149 :type n_max_naive: :class:`int`
150 :param oum: Override USalign Mapping. Inject rigid_mapping of
151 :class:`Scorer` object into USalign to compute TM-score.
152 Experimental feature with limitations.
153 :type oum: :class:`bool`
154 :param min_pep_length: Relevant parameter if short peptides are involved in
155 scoring. Minimum peptide length for a chain in the
156 target structure to be considered in chain mapping.
157 The chain mapping algorithm first performs an all vs.
158 all pairwise sequence alignment to identify \"equal\"
159 chains within the target structure. We go for simple
160 sequence identity there. Short sequences can be
161 problematic as they may produce high sequence identity
162 alignments by pure chance.
163 :type min_pep_length: :class:`int`
164 :param min_nuc_length: Relevant parameter if short nucleotides are involved
165 in scoring. Minimum nucleotide length for a chain in
166 the target structure to be considered in chain
167 mapping. The chain mapping algorithm first performs
168 an all vs. all pairwise sequence alignment to
169 identify \"equal\" chains within the target
170 structure. We go for simple sequence identity there.
171 Short sequences can be problematic as they may
172 produce high sequence identity alignments by pure
174 :type min_nuc_length: :class:`int`
175 :param lddt_add_mdl_contacts: lDDT specific flag. Only using contacts in
176 lDDT that are within a certain distance
177 threshold in the target does not penalize
178 for added model contacts. If set to True, this
179 flag will also consider target contacts
180 that are within the specified distance
181 threshold in the model but not necessarily in
182 the target. No contact will be added if the
183 respective atom pair is not resolved in the
185 :type lddt_add_mdl_contacts: :class:`bool`
186 :param dockq_capri_peptide: Flag that changes two things in the way DockQ
187 and its underlying scores are computed which is
188 proposed by the CAPRI community when scoring
189 peptides (PMID: 31886916).
190 ONE: Two residues are considered in contact if
191 any of their atoms is within 5A. This is
192 relevant for fnat and fnonat scores. CAPRI
193 suggests to lower this threshold to 4A for
194 protein-peptide interactions.
195 TWO: irmsd is computed on interface residues. A
196 residue is defined as interface residue if any
197 of its atoms is within 10A of another chain.
198 CAPRI suggests to lower the default of 10A to
199 8A in combination with only considering CB atoms
200 for protein-peptide interactions.
201 This flag has no influence on patch_dockq
203 :type dockq_capri_peptide: :class:`bool`
204 :param lddt_symmetry_settings: Passed as *symmetry_settings* parameter to
205 lDDT scorer. Default: None
206 :type lddt_symmetry_settings: :class:`ost.mol.alg.lddt.SymmetrySettings`
207 :param lddt_inclusion_radius: lDDT inclusion radius.
209 def __init__(self, model, target, resnum_alignments=False,
210 molck_settings = None, cad_score_exec = None,
211 custom_mapping=None, custom_rigid_mapping=None,
212 usalign_exec = None, lddt_no_stereochecks=False,
213 n_max_naive=40320, oum=False, min_pep_length = 6,
214 min_nuc_length = 4, lddt_add_mdl_contacts=False,
215 dockq_capri_peptide=False, lddt_symmetry_settings = None,
216 lddt_inclusion_radius = 15.0):
231 if molck_settings
is None:
236 rm_zero_occ_atoms=
False,
240 LogScript(
"Cleaning up input structures")
241 Molck(self.
_model_model, conop.GetDefaultLib(), molck_settings)
242 Molck(self.
_target_target, conop.GetDefaultLib(), molck_settings)
244 if resnum_alignments:
257 self.
_model_model = self.
_model_model.Select(
"peptide=True or nucleotide=True")
258 self.
_target_target = self.
_target_target.Select(
"peptide=True or nucleotide=True")
261 for ch
in self.
_model_model.chains:
262 if ch.GetName().strip() ==
"":
263 raise RuntimeError(
"Model chains must have valid chain names")
264 if ch.GetName().strip() ==
"'" or ch.GetName().strip() ==
'"':
265 raise RuntimeError(
"Model chains cannot be named with quote signs (' or \"\")")
268 for ch
in self.
_target_target.chains:
269 if ch.GetName().strip() ==
"":
270 raise RuntimeError(
"Target chains must have valid chain names")
271 if ch.GetName().strip() ==
"'" or ch.GetName().strip() ==
'"':
272 raise RuntimeError(
"Target chains cannot be named with quote signs (' or \"\")")
274 if resnum_alignments:
278 for ch
in self.
_model_model.chains:
279 ins_codes = [r.GetNumber().GetInsCode()
for r
in ch.residues]
280 if len(set(ins_codes)) != 1
or ins_codes[0] !=
'\0':
281 raise RuntimeError(
"Residue numbers in each model chain "
282 "must not contain insertion codes if "
283 "resnum_alignments are enabled")
284 nums = [r.GetNumber().GetNum()
for r
in ch.residues]
285 if not all(i < j
for i, j
in zip(nums, nums[1:])):
286 raise RuntimeError(
"Residue numbers in each model chain "
287 "must be strictly increasing if "
288 "resnum_alignments are enabled")
290 for ch
in self.
_target_target.chains:
291 ins_codes = [r.GetNumber().GetInsCode()
for r
in ch.residues]
292 if len(set(ins_codes)) != 1
or ins_codes[0] !=
'\0':
293 raise RuntimeError(
"Residue numbers in each target chain "
294 "must not contain insertion codes if "
295 "resnum_alignments are enabled")
296 nums = [r.GetNumber().GetNum()
for r
in ch.residues]
297 if not all(i < j
for i, j
in zip(nums, nums[1:])):
298 raise RuntimeError(
"Residue numbers in each target chain "
299 "must be strictly increasing if "
300 "resnum_alignments are enabled")
302 if usalign_exec
is not None:
303 if not os.path.exists(usalign_exec):
304 raise RuntimeError(f
"USalign exec ({usalign_exec}) "
306 if not os.access(usalign_exec, os.X_OK):
307 raise RuntimeError(f
"USalign exec ({usalign_exec}) "
308 f
"is not executable")
351 self.
_lddt_lddt =
None
402 self.
_fnat_fnat =
None
406 self.
_nnat_nnat =
None
407 self.
_nmdl_nmdl =
None
438 self.
_rmsd_rmsd =
None
449 if custom_mapping
is not None:
452 if custom_rigid_mapping
is not None:
455 LogDebug(
"Scorer sucessfully initialized")
459 """ Model with Molck cleanup
461 :type: :class:`ost.mol.EntityHandle`
467 """ The original model passed at object construction
469 :type: :class:`ost.mol.EntityHandle`/:class:`ost.mol.EntityView`
475 """ Target with Molck cleanup
477 :type: :class:`ost.mol.EntityHandle`
483 """ The original target passed at object construction
485 :type: :class:`ost.mol.EntityHandle`/:class:`ost.mol.EntityView`
491 """ Alignments of :attr:`~model`/:attr:`~target` chains
493 Alignments for each pair of chains mapped in :attr:`~mapping`.
494 First sequence is target sequence, second sequence the model sequence.
496 :type: :class:`list` of :class:`ost.seq.AlignmentHandle`
498 if self.
_aln_aln
is None:
504 """ Stereochecked equivalent of :attr:`~aln`
506 The alignments may differ, as stereochecks potentially remove residues
508 :type: :class:`list` of :class:`ost.seq.AlignmentHandle`
516 """ Alignments of :attr:`~model_orig`/:attr:`~target_orig` chains
518 Selects for peptide and nucleotide residues before sequence
519 extraction. Includes residues that would be removed by molck in
520 structure preprocessing.
522 :type: :class:`list` of :class:`ost.seq.AlignmentHandle`
530 """ Alignments of :attr:`~trimmed_model`/:attr:`~target` chains
532 Alignments for each pair of chains mapped in :attr:`~mapping`.
533 First sequence is target sequence, second sequence the model sequence.
535 :type: :class:`list` of :class:`ost.seq.AlignmentHandle`
543 """ View of :attr:`~model` that has stereochemistry checks applied
545 First, a selection for peptide/nucleotide residues is performed,
546 secondly peptide sidechains with stereochemical irregularities are
547 removed (full residue if backbone atoms are involved). Irregularities
548 are clashes or bond lengths/angles more than 12 standard deviations
549 from expected values.
551 :type: :class:`ost.mol.EntityView`
559 """ Clashing model atoms
561 :type: :class:`list` of :class:`ost.mol.alg.stereochemistry.ClashInfo`
569 """ Model bonds with unexpected stereochemistry
571 :type: :class:`list` of
572 :class:`ost.mol.alg.stereochemistry.BondViolationInfo`
580 """ Model angles with unexpected stereochemistry
582 :type: :class:`list` of
583 :class:`ost.mol.alg.stereochemistry.AngleViolationInfo`
591 """ Same as :attr:`~stereochecked_model` for :attr:`~target`
593 :type: :class:`ost.mol.EntityView`
601 """ Clashing target atoms
603 :type: :class:`list` of :class:`ost.mol.alg.stereochemistry.ClashInfo`
611 """ Target bonds with unexpected stereochemistry
613 :type: :class:`list` of
614 :class:`ost.mol.alg.stereochemistry.BondViolationInfo`
622 """ Target angles with unexpected stereochemistry
624 :type: :class:`list` of
625 :class:`ost.mol.alg.stereochemistry.AngleViolationInfo`
633 """ :attr:`~model` trimmed to target
635 Removes residues that are not covered by :class:`target` given
636 :attr:`~mapping`. In other words: no model residues without experimental
637 evidence from :class:`target`.
639 :type: :class:`ost.mol.EntityView`
647 """ Chain mapper object for given :attr:`~target`
649 :type: :class:`ost.mol.alg.chain_mapping.ChainMapper`
661 """ Full chain mapping result for :attr:`~target`/:attr:`~model`
663 Computed with :func:`ost.mol.alg.ChainMapper.GetMapping`
665 :type: :class:`ost.mol.alg.chain_mapping.MappingResult`
668 LogScript(
"Computing chain mapping")
676 """ Full chain mapping result for :attr:`~target`/:attr:`~model`
678 Computed with :func:`ost.mol.alg.ChainMapper.GetRMSDMapping`
680 :type: :class:`ost.mol.alg.chain_mapping.MappingResult`
683 LogScript(
"Computing rigid chain mapping")
690 """ Interface residues in :attr:`~model`
692 Thats all residues having a contact with at least one residue from
693 another chain (CB-CB distance <= 8A, CA in case of Glycine)
695 :type: :class:`dict` with chain names as key and and :class:`list`
696 with residue numbers of the respective interface residues.
705 """ Same as :attr:`~model_interface_residues` for :attr:`~target`
707 :type: :class:`dict` with chain names as key and and :class:`list`
708 with residue numbers of the respective interface residues.
717 """ lDDT scorer for :attr:`~target`/:attr:`~stereochecked_target`
719 Depending on :attr:`~lddt_no_stereocheck` and
720 :attr:`~lddt_symmetry_settings`.
722 :type: :class:`ost.mol.alg.lddt.lDDTScorer`
737 """ Backbone only lDDT scorer for :attr:`~target`
739 No stereochecks applied for bb only lDDT which considers CA atoms
740 for peptides and C3' atoms for nucleotides.
742 :type: :class:`ost.mol.alg.lddt.lDDTScorer`
752 """ QS scorer constructed from :attr:`~mapping`
754 The scorer object is constructed with default parameters and relates to
755 :attr:`~model` and :attr:`~target` (no stereochecks).
757 :type: :class:`ost.mol.alg.qsscore.QSScorer`
773 self.
mappingmapping.chem_groups,
780 """ Global lDDT score in range [0.0, 1.0]
782 Computed based on :attr:`~stereochecked_model`. In case of oligomers,
783 :attr:`~mapping` is used.
785 :type: :class:`float`
787 if self.
_lddt_lddt
is None:
789 return self.
_lddt_lddt
793 """ Per residue lDDT scores in range [0.0, 1.0]
795 Computed based on :attr:`~stereochecked_model` but scores for all
796 residues in :attr:`~model` are reported. If a residue has been removed
797 by stereochemistry checks, the respective score is set to 0.0. If a
798 residue is not covered by the target or is in a chain skipped by the
799 chain mapping procedure (happens for super short chains), the respective
800 score is set to None. In case of oligomers, :attr:`~mapping` is used.
810 """ Per atom lDDT scores in range [0.0, 1.0]
812 Computed based on :attr:`~stereochecked_model` but scores for all
813 atoms in :attr:`~model` are reported. If an atom has been removed
814 by stereochemistry checks, the respective score is set to 0.0. If an
815 atom is not covered by the target or is in a chain skipped by the
816 chain mapping procedure (happens for super short chains), the respective
817 score is set to None. In case of oligomers, :attr:`~mapping` is used.
827 """ Backbone only global lDDT score in range [0.0, 1.0]
829 Computed based on :attr:`~model` on backbone atoms only. This is CA for
830 peptides and C3' for nucleotides. No stereochecks are performed. In case
831 of oligomers, :attr:`~mapping` is used.
833 :type: :class:`float`
841 """ Backbone only per residue lDDT scores in range [0.0, 1.0]
843 Computed based on :attr:`~model` on backbone atoms only. This is CA for
844 peptides and C3' for nucleotides. No stereochecks are performed. If a
845 residue is not covered by the target or is in a chain skipped by the
846 chain mapping procedure (happens for super short chains), the respective
847 score is set to None. In case of oligomers, :attr:`~mapping` is used.
857 """ Global interface lDDT score in range [0.0, 1.0]
859 This is lDDT only based on inter-chain contacts. Value is None if no
860 such contacts are present. For example if we're dealing with a monomer.
861 Computed based on :attr:`~stereochecked_model` and :attr:`~mapping` for
864 :type: :class:`float`
866 if self.
_ilddt_ilddt
is None:
877 Computed based on :attr:`~model` using :attr:`~mapping`
879 :type: :class:`float`
887 """ Global QS-score - only computed on aligned residues
889 Computed based on :attr:`~model` using :attr:`~mapping`. The QS-score
890 computation only considers contacts between residues with a mapping
891 between target and model. As a result, the score won't be lowered in
892 case of additional chains/residues in any of the structures.
894 :type: :class:`float`
902 """ Interfaces in :attr:`~target` with non-zero contribution to
903 :attr:`~qs_global`/:attr:`~qs_best`
905 Chain names are lexicographically sorted.
907 :type: :class:`list` of :class:`tuple` with 2 elements each:
918 """ Interfaces in :attr:`~model` with non-zero contribution to
919 :attr:`~qs_global`/:attr:`~qs_best`
921 Chain names are lexicographically sorted.
923 :type: :class:`list` of :class:`tuple` with 2 elements each:
935 """ Interfaces in :attr:`~qs_target_interfaces` that can be mapped
938 Target chain names are lexicographically sorted.
940 :type: :class:`list` of :class:`tuple` with 4 elements each:
941 (trg_ch1, trg_ch2, mdl_ch1, mdl_ch2)
945 flat_mapping = self.
mappingmapping.GetFlatMapping()
947 if i[0]
in flat_mapping
and i[1]
in flat_mapping:
955 """ QS-score for each interface in :attr:`~qs_interfaces`
957 :type: :class:`list` of :class:`float`
965 """ QS-score for each interface in :attr:`~qs_interfaces`
967 Only computed on aligned residues
969 :type: :class:`list` of :class:`float`
979 A contact is a pair or residues from distinct chains that have
980 a minimal heavy atom distance < 5A. Contacts are specified as
981 :class:`tuple` with two strings in format:
982 <cname>.<rnum>.<ins_code>
984 :type: :class:`list` of :class:`tuple`
992 """ Same for :attr:`~model`
1000 """ Same for :attr:`~trimmed_model`
1008 """ Interfaces in :class:`target` which have at least one contact
1010 Contact as defined in :attr:`~native_contacts`,
1011 chain names are lexicographically sorted.
1013 :type: :class:`list` of :class:`tuple` with 2 elements each
1018 tmp = [(min(x[0],x[1]), max(x[0],x[1]))
for x
in tmp]
1024 """ Interfaces in :class:`model` which have at least one contact
1026 Contact as defined in :attr:`~native_contacts`,
1027 chain names are lexicographically sorted.
1029 :type: :class:`list` of :class:`tuple` with 2 elements each
1034 tmp = [(min(x[0],x[1]), max(x[0],x[1]))
for x
in tmp]
1040 """ Fraction of model contacts that are also present in target
1042 :type: :class:`float`
1050 """ Fraction of target contacts that are correctly reproduced in model
1052 :type: :class:`float`
1060 """ ICS (Interface Contact Similarity) score
1062 Combination of :attr:`~ics_precision` and :attr:`~ics_recall`
1063 using the F1-measure
1065 :type: :class:`float`
1067 if self.
_ics_ics
is None:
1069 return self.
_ics_ics
1073 """ Per-interface ICS precision
1075 :attr:`~ics_precision` for each interface in
1076 :attr:`~contact_target_interfaces`
1078 :type: :class:`list` of :class:`float`
1087 """ Per-interface ICS recall
1089 :attr:`~ics_recall` for each interface in
1090 :attr:`~contact_target_interfaces`
1092 :type: :class:`list` of :class:`float`
1100 """ Per-interface ICS (Interface Contact Similarity) score
1102 :attr:`~ics` for each interface in
1103 :attr:`~contact_target_interfaces`
1105 :type: :class:`float`
1114 """ Fraction of model interface residues that are also interface
1117 :type: :class:`float`
1125 """ Fraction of target interface residues that are also interface
1128 :type: :class:`float`
1136 """ IPS (Interface Patch Similarity) score
1138 Jaccard coefficient of interface residues in target and their mapped
1139 counterparts in model
1141 :type: :class:`float`
1143 if self.
_ips_ips
is None:
1145 return self.
_ips_ips
1149 """ Same as :attr:`~ics` but with trimmed model
1151 Model is trimmed to residues which can me mapped to target in order
1152 to not penalize contacts in the model for which we have no experimental
1155 :type: :class:`float`
1163 """ Same as :attr:`~ics_precision` but with trimmed model
1165 Model is trimmed to residues which can me mapped to target in order
1166 to not penalize contacts in the model for which we have no experimental
1169 :type: :class:`float`
1177 """ Same as :attr:`~ics_recall` but with trimmed model
1179 Model is trimmed to residues which can me mapped to target in order
1180 to not penalize contacts in the model for which we have no experimental
1183 :type: :class:`float`
1191 """ Same as :attr:`~per_interface_ics_precision` but with :attr:`~trimmed_model`
1193 :attr:`~ics_precision_trimmed` for each interface in
1194 :attr:`~contact_target_interfaces`
1196 :type: :class:`list` of :class:`float`
1205 """ Same as :attr:`~per_interface_ics_recall` but with :attr:`~trimmed_model`
1207 :attr:`~ics_recall_trimmed` for each interface in
1208 :attr:`~contact_target_interfaces`
1210 :type: :class:`list` of :class:`float`
1218 """ Same as :attr:`~per_interface_ics` but with :attr:`~trimmed_model`
1220 :attr:`~ics` for each interface in
1221 :attr:`~contact_target_interfaces`
1223 :type: :class:`float`
1232 """ Same as :attr:`~ips` but with trimmed model
1234 Model is trimmed to residues which can me mapped to target in order
1235 to not penalize contacts in the model for which we have no experimental
1238 :type: :class:`float`
1246 """ Same as :attr:`~ips_precision` but with trimmed model
1248 Model is trimmed to residues which can me mapped to target in order
1249 to not penalize contacts in the model for which we have no experimental
1252 :type: :class:`float`
1260 """ Same as :attr:`~ips_recall` but with trimmed model
1262 Model is trimmed to residues which can me mapped to target in order
1263 to not penalize contacts in the model for which we have no experimental
1266 :type: :class:`float`
1274 """ Per-interface IPS precision
1276 :attr:`~ips_precision` for each interface in
1277 :attr:`~contact_target_interfaces`
1279 :type: :class:`list` of :class:`float`
1287 """ Per-interface IPS recall
1289 :attr:`~ips_recall` for each interface in
1290 :attr:`~contact_target_interfaces`
1292 :type: :class:`list` of :class:`float`
1300 """ Per-interface IPS (Interface Patch Similarity) score
1302 :attr:`~ips` for each interface in
1303 :attr:`~contact_target_interfaces`
1305 :type: :class:`list` of :class:`float`
1314 """ Same as :attr:`~per_interface_ips_precision` but with :attr:`~trimmed_model`
1316 :attr:`~ips_precision_trimmed` for each interface in
1317 :attr:`~contact_target_interfaces`
1319 :type: :class:`list` of :class:`float`
1328 """ Same as :attr:`~per_interface_ips_recall` but with :attr:`~trimmed_model`
1330 :attr:`~ics_recall_trimmed` for each interface in
1331 :attr:`~contact_target_interfaces`
1333 :type: :class:`list` of :class:`float`
1341 """ Same as :attr:`~per_interface_ips` but with :attr:`~trimmed_model`
1343 :attr:`~ics` for each interface in
1344 :attr:`~contact_target_interfaces`
1346 :type: :class:`float`
1355 """ Interfaces in :attr:`~target` that are relevant for DockQ
1357 All interfaces in :attr:`~target` with non-zero contacts that are
1358 relevant for DockQ. Includes protein-protein, protein-nucleotide and
1359 nucleotide-nucleotide interfaces. Chain names for each interface are
1360 lexicographically sorted.
1362 :type: :class:`list` of :class:`tuple` with 2 elements each:
1372 contact_d = contact_d)
1375 interfaces = cent.interacting_chains
1376 interfaces = [(min(x[0],x[1]), max(x[0],x[1]))
for x
in interfaces]
1378 pep_seqs = set([s.GetName()
for s
in self.
chain_mapperchain_mapper.polypep_seqs])
1379 nuc_seqs = set([s.GetName()
for s
in self.
chain_mapperchain_mapper.polynuc_seqs])
1381 seqs = pep_seqs.union(nuc_seqs)
1383 for interface
in interfaces:
1384 if interface[0]
in seqs
and interface[1]
in seqs:
1391 """ Interfaces in :attr:`~dockq_target_interfaces` that can be mapped
1394 Target chain names are lexicographically sorted
1396 :type: :class:`list` of :class:`tuple` with 4 elements each:
1397 (trg_ch1, trg_ch2, mdl_ch1, mdl_ch2)
1401 flat_mapping = self.
mappingmapping.GetFlatMapping()
1403 if i[0]
in flat_mapping
and i[1]
in flat_mapping:
1406 flat_mapping[i[1]]))
1411 """ DockQ scores for interfaces in :attr:`~dockq_interfaces`
1413 :class:`list` of :class:`float`
1421 """ fnat scores for interfaces in :attr:`~dockq_interfaces`
1423 fnat: Fraction of native contacts that are also present in model
1425 :class:`list` of :class:`float`
1427 if self.
_fnat_fnat
is None:
1429 return self.
_fnat_fnat
1433 """ N native contacts for interfaces in :attr:`~dockq_interfaces`
1435 :class:`list` of :class:`int`
1437 if self.
_nnat_nnat
is None:
1439 return self.
_nnat_nnat
1443 """ N model contacts for interfaces in :attr:`~dockq_interfaces`
1445 :class:`list` of :class:`int`
1447 if self.
_nmdl_nmdl
is None:
1449 return self.
_nmdl_nmdl
1453 """ fnonnat scores for interfaces in :attr:`~dockq_interfaces`
1455 fnat: Fraction of model contacts that are not present in target
1457 :class:`list` of :class:`float`
1465 """ irmsd scores for interfaces in :attr:`~dockq_interfaces`
1467 irmsd: RMSD of interface (RMSD computed on backbone atoms) which
1468 consists of each residue that has at least one heavy atom within 10A of
1469 other chain. Backbone atoms for proteins: "CA","C","N","O", for
1470 nucleotides: "P", "OP1", "OP2", "O2'", "O3'", "O4'", "O5'", "C1'",
1471 "C2'", "C3'", "C4'", "C5'".
1473 :class:`list` of :class:`float`
1475 if self.
_irmsd_irmsd
is None:
1481 """ lrmsd scores for interfaces in :attr:`~dockq_interfaces`
1483 lrmsd: The two chains involved in the interface are superposed based on
1484 the receptor (rigid min RMSD superposition) and the ligand RMSD is
1485 reported. Receptor is the chain with more residues. Superposition and
1486 RMSD is computed on same backbone atoms as :attr:`~irmsd`.
1488 :class:`list` of :class:`float`
1490 if self.
_lrmsd_lrmsd
is None:
1496 """ Average of DockQ scores in :attr:`~dockq_scores`
1498 In its original implementation, DockQ only operates on single
1499 interfaces. Thus the requirement to combine scores for higher order
1502 :type: :class:`float`
1510 """ Same as :attr:`~dockq_ave`, weighted by native contacts
1512 :type: :class:`float`
1520 """ Same as :attr:`~dockq_ave` but penalizing for missing interfaces
1522 Interfaces that are not covered in model are added as 0.0
1523 in average computation.
1525 :type: :class:`float`
1533 """ Same as :attr:`~dockq_ave_full`, but weighted
1535 Interfaces that are not covered in model are added as 0.0 in
1536 average computations and the respective weights are derived from
1537 number of contacts in respective target interface.
1545 """ Mapped representative positions in target
1547 Thats CA positions for peptide residues and C3' positions for
1548 nucleotides. Has same length as :attr:`~mapped_model_pos` and mapping
1549 is based on :attr:`~mapping`.
1551 :type: :class:`ost.geom.Vec3List`
1559 """ Mapped representative positions in target
1561 Thats the equivalent of :attr:`~mapped_target_pos` but containing more
1562 backbone atoms (N, CA, C for peptide residues and O5', C5', C4', C3', O3
1563 for nucleotide residues). mapping is based on :attr:`~mapping`.
1565 :type: :class:`ost.geom.Vec3List`
1573 """ Mapped representative positions in model
1575 Thats CA positions for peptide residues and C3' positions for
1576 nucleotides. Has same length as :attr:`~mapped_target_pos` and mapping
1577 is based on :attr:`~mapping`.
1579 :type: :class:`ost.geom.Vec3List`
1587 """ Mapped representative positions in model
1589 Thats the equivalent of :attr:`~mapped_model_pos` but containing more
1590 backbone atoms (N, CA, C for peptide residues and O5', C5', C4', C3', O3
1591 for nucleotide residues). mapping is based on :attr:`~mapping`.
1593 :type: :class:`ost.geom.Vec3List`
1601 """ :attr:`~mapped_model_pos` with :attr:`~transform` applied
1603 :type: :class:`ost.geom.Vec3List`
1613 """ Number of target residues which have no mapping to model
1623 """ Transform: :attr:`~mapped_model_pos` onto :attr:`~mapped_target_pos`
1625 Computed using Kabsch minimal rmsd algorithm. If number of positions
1626 is too small (< 3), :attr:`~mapped_model_pos_full_bb` and
1627 :attr:`~mapped_target_pos_full_bb` are used.
1629 :type: :class:`ost.geom.Mat4`
1636 self.
_transform_transform = res.transformation
1643 self.
_transform_transform = res.transformation
1648 """ Mapped representative positions in target
1650 Thats CA positions for peptide residues and C3' positions for
1651 nucleotides. Has same length as :attr:`~rigid_mapped_model_pos` and mapping
1652 is based on :attr:`~rigid_mapping`.
1654 :type: :class:`ost.geom.Vec3List`
1662 """ Mapped representative positions in target
1664 Thats the equivalent of :attr:`~rigid_mapped_target_pos` but containing
1665 more backbone atoms (N, CA, C for peptide residues and O5', C5', C4',
1666 C3', O3 for nucleotide residues). mapping is based on :attr:`~mapping`.
1668 :type: :class:`ost.geom.Vec3List`
1676 """ Mapped representative positions in model
1678 Thats CA positions for peptide residues and C3' positions for
1679 nucleotides. Has same length as :attr:`~mapped_target_pos` and mapping
1680 is based on :attr:`~rigid_mapping`.
1682 :type: :class:`ost.geom.Vec3List`
1690 """ Mapped representative positions in model
1692 Thats the equivalent of :attr:`~rigid_mapped_model_pos` but containing
1693 more backbone atoms (N, CA, C for peptide residues and O5', C5', C4',
1694 C3', O3 for nucleotide residues). mapping is based on :attr:`~mapping`.
1696 :type: :class:`ost.geom.Vec3List`
1704 """ :attr:`~rigid_mapped_model_pos` with :attr:`~rigid_transform` applied
1706 :type: :class:`ost.geom.Vec3List`
1716 """ Number of target residues which have no rigid mapping to model
1726 """ Transform: :attr:`~rigid_mapped_model_pos` onto :attr:`~rigid_mapped_target_pos`
1728 Computed using Kabsch minimal rmsd algorithm. If number of positions
1729 is too small (< 3), :attr:`~rigid_mapped_model_pos_full_bb` and
1730 :attr:`~rigid_mapped_target_pos_full_bb` are used.
1732 :type: :class:`ost.geom.Mat4`
1751 """ Fraction CA (C3' for nucleotides) that can be superposed within 0.5A
1753 Uses :attr:`~rigid_mapped_model_pos` and :attr:`~rigid_mapped_target_pos`.
1754 Similar iterative algorithm as LGA tool
1756 :type: :class:`float`
1758 if self.
_gdt_05_gdt_05
is None:
1761 for window_size
in wsizes:
1764 window_size, 1000, 0.5)[0]
1769 self.
_gdt_05_gdt_05 = float(n) / n_full
1776 """ Fraction CA (C3' for nucleotides) that can be superposed within 1.0A
1778 Uses :attr:`~rigid_mapped_model_pos` and :attr:`~rigid_mapped_target_pos`.
1779 Similar iterative algorithm as LGA tool
1781 :type: :class:`float`
1783 if self.
_gdt_1_gdt_1
is None:
1786 for window_size
in wsizes:
1789 window_size, 1000, 1.0)[0]
1794 self.
_gdt_1_gdt_1 = float(n) / n_full
1801 """ Fraction CA (C3' for nucleotides) that can be superposed within 2.0A
1803 Uses :attr:`~rigid_mapped_model_pos` and :attr:`~rigid_mapped_target_pos`.
1804 Similar iterative algorithm as LGA tool
1807 :type: :class:`float`
1809 if self.
_gdt_2_gdt_2
is None:
1812 for window_size
in wsizes:
1815 window_size, 1000, 2.0)[0]
1820 self.
_gdt_2_gdt_2 = float(n) / n_full
1827 """ Fraction CA (C3' for nucleotides) that can be superposed within 4.0A
1829 Uses :attr:`~rigid_mapped_model_pos` and :attr:`~rigid_mapped_target_pos`.
1830 Similar iterative algorithm as LGA tool
1832 :type: :class:`float`
1834 if self.
_gdt_4_gdt_4
is None:
1837 for window_size
in wsizes:
1840 window_size, 1000, 4.0)[0]
1845 self.
_gdt_4_gdt_4 = float(n) / n_full
1852 """ Fraction CA (C3' for nucleotides) that can be superposed within 8.0A
1854 Similar iterative algorithm as LGA tool
1856 :type: :class:`float`
1858 if self.
_gdt_8_gdt_8
is None:
1861 for window_size
in wsizes:
1864 window_size, 1000, 8.0)[0]
1869 self.
_gdt_8_gdt_8 = float(n) / n_full
1877 """ avg GDT with thresholds: 8.0A, 4.0A, 2.0A and 1.0A
1879 :type: :class:`float`
1881 if self.
_gdtts_gdtts
is None:
1882 LogScript(
"Computing GDT-TS score")
1888 """ avg GDT with thresholds: 4.0A, 2.0A, 1.0A and 0.5A
1890 :type: :class:`float`
1892 if self.
_gdtha_gdtha
is None:
1893 LogScript(
"Computing GDT-HA score")
1901 Computed on :attr:`~rigid_transformed_mapped_model_pos` and
1902 :attr:`~rigid_mapped_target_pos`
1904 :type: :class:`float`
1906 if self.
_rmsd_rmsd
is None:
1907 LogScript(
"Computing RMSD")
1910 return self.
_rmsd_rmsd
1914 """ The global CAD atom-atom (AA) score
1916 Computed based on :attr:`~model`. In case of oligomers, :attr:`~mapping`
1919 :type: :class:`float`
1927 """ The per-residue CAD atom-atom (AA) scores
1929 Computed based on :attr:`~model`. In case of oligomers, :attr:`~mapping`
1932 :type: :class:`dict`
1940 """ Patch QS-scores for each residue in :attr:`~model_interface_residues`
1942 Representative patches for each residue r in chain c are computed as
1945 * mdl_patch_one: All residues in c with CB (CA for GLY) positions within
1946 8A of r and within 12A of residues from any other chain.
1947 * mdl_patch_two: Closest residue x to r in any other chain gets
1948 identified. Patch is then constructed by selecting all residues from
1949 any other chain within 8A of x and within 12A from any residue in c.
1950 * trg_patch_one: Chain name and residue number based mapping from
1952 * trg_patch_two: Chain name and residue number based mapping from
1955 Results are stored in the same manner as
1956 :attr:`~model_interface_residues`, with corresponding scores instead of
1957 residue numbers. Scores for residues which are not
1958 :class:`mol.ChemType.AMINOACIDS` are set to None. Additionally,
1959 interface patches are derived from :attr:`~model`. If they contain
1960 residues which are not covered by :attr:`~target`, the score is set to
1963 :type: :class:`dict` with chain names as key and and :class:`list`
1964 with scores of the respective interface residues.
1972 """ Same as :attr:`~patch_qs` but for DockQ scores
1980 """ TM-score computed with USalign
1982 USalign executable can be specified with usalign_exec kwarg at Scorer
1983 construction, an OpenStructure internal copy of the USalign code is
1986 :type: :class:`float`
1994 """ Mapping computed with USalign
1996 Dictionary with target chain names as key and model chain names as
1997 values. No guarantee that all chains are mapped. USalign executable
1998 can be specified with usalign_exec kwarg at Scorer construction, an
1999 OpenStructure internal copy of the USalign code is used otherwise.
2001 :type: :class:`dict`
2007 def _aln_helper(self, target, model):
2012 for ch
in target.chains:
2013 cname = ch.GetName()
2014 s =
''.join([r.one_letter_code
for r
in ch.residues])
2015 s = seq.CreateSequence(ch.GetName(), s)
2016 s.AttachView(target.Select(f
"cname={mol.QueryQuoteName(cname)}"))
2017 trg_seqs[ch.GetName()] = s
2019 for ch
in model.chains:
2020 cname = ch.GetName()
2021 s =
''.join([r.one_letter_code
for r
in ch.residues])
2022 s = seq.CreateSequence(cname, s)
2023 s.AttachView(model.Select(f
"cname={mol.QueryQuoteName(cname)}"))
2024 mdl_seqs[ch.GetName()] = s
2027 trg_pep_chains = [s.GetName()
for s
in self.
chain_mapperchain_mapper.polypep_seqs]
2028 trg_nuc_chains = [s.GetName()
for s
in self.
chain_mapperchain_mapper.polynuc_seqs]
2029 trg_pep_chains = set(trg_pep_chains)
2030 trg_nuc_chains = set(trg_nuc_chains)
2031 for trg_ch, mdl_ch
in self.
mappingmapping.GetFlatMapping().items():
2032 if mdl_ch
in mdl_seqs
and trg_ch
in trg_seqs:
2033 if trg_ch
in trg_pep_chains:
2034 stype = mol.ChemType.AMINOACIDS
2035 elif trg_ch
in trg_nuc_chains:
2036 stype = mol.ChemType.NUCLEOTIDES
2038 raise RuntimeError(
"Chain name inconsistency... ask "
2040 alns.append(self.
chain_mapperchain_mapper.Align(trg_seqs[trg_ch],
2043 alns[-1].AttachView(0, trg_seqs[trg_ch].GetAttachedView())
2044 alns[-1].AttachView(1, mdl_seqs[mdl_ch].GetAttachedView())
2047 def _compute_pepnuc_aln(self):
2048 query =
"peptide=true or nucleotide=true"
2049 pep_nuc_target = self.
target_origtarget_orig.Select(query)
2050 pep_nuc_model = self.
model_origmodel_orig.Select(query)
2053 def _compute_aln(self):
2056 def _compute_stereochecked_aln(self):
2060 def _compute_lddt(self):
2061 LogScript(
"Computing all-atom lDDT")
2063 flat_mapping = self.
mappingmapping.GetFlatMapping(mdl_as_key=
True)
2066 stereochecked_alns = dict()
2068 mdl_seq = aln.GetSequence(1)
2069 stereochecked_alns[mdl_seq.name] = aln
2071 for aln
in self.
alnaln:
2072 mdl_seq = aln.GetSequence(1)
2073 alns[mdl_seq.name] = aln
2078 aa_local_lddt =
None
2081 lddt_chain_mapping = dict()
2082 for mdl_ch, trg_ch
in flat_mapping.items():
2084 lddt_chain_mapping[mdl_ch] = trg_ch
2086 chain_mapping = lddt_chain_mapping,
2087 residue_mapping = alns,
2088 check_resnames=
False,
2089 local_lddt_prop=
"lddt",
2091 set_atom_props=
True)[0]
2093 aa_local_lddt = dict()
2094 for r
in self.
modelmodel.residues:
2096 cname = r.GetChain().GetName()
2097 if cname
not in local_lddt:
2098 local_lddt[cname] = dict()
2099 aa_local_lddt[cname] = dict()
2101 rnum = r.GetNumber()
2102 if rnum
not in aa_local_lddt[cname]:
2103 aa_local_lddt[cname][rnum] = dict()
2105 if r.HasProp(
"lddt"):
2106 score = round(r.GetFloatProp(
"lddt"), 3)
2107 local_lddt[cname][rnum] = score
2111 local_lddt[cname][rnum] =
None
2114 if a.HasProp(
"lddt"):
2115 score = round(a.GetFloatProp(
"lddt"), 3)
2116 aa_local_lddt[cname][rnum][a.GetName()] = score
2121 aa_local_lddt[cname][rnum][a.GetName()] =
None
2125 lddt_chain_mapping = dict()
2126 for mdl_ch, trg_ch
in flat_mapping.items():
2127 if mdl_ch
in stereochecked_alns:
2128 lddt_chain_mapping[mdl_ch] = trg_ch
2130 chain_mapping = lddt_chain_mapping,
2131 residue_mapping = stereochecked_alns,
2132 check_resnames=
False,
2133 local_lddt_prop=
"lddt",
2135 set_atom_props=
True)[0]
2137 aa_local_lddt = dict()
2138 for r
in self.
modelmodel.residues:
2139 cname = r.GetChain().GetName()
2140 if cname
not in local_lddt:
2141 local_lddt[cname] = dict()
2142 aa_local_lddt[cname] = dict()
2143 rnum = r.GetNumber()
2144 if rnum
not in aa_local_lddt[cname]:
2145 aa_local_lddt[cname][rnum] = dict()
2147 if r.HasProp(
"lddt"):
2148 score = round(r.GetFloatProp(
"lddt"), 3)
2149 local_lddt[cname][rnum] = score
2155 if a.HasProp(
"lddt"):
2156 score = round(a.GetFloatProp(
"lddt"), 3)
2157 aa_local_lddt[cname][rnum][a.GetName()] = score
2167 if cname
in flat_mapping:
2168 for col
in alns[cname]:
2169 if col[0] !=
'-' and col[1] !=
'-':
2170 if col.GetResidue(1).number == r.number:
2171 trg_r = col.GetResidue(0)
2173 if trg_r
is not None:
2174 trg_cname = trg_r.GetChain().GetName()
2175 trg_rnum = trg_r.GetNumber()
2186 if trg_r
is not None and not trg_r.FindAtom(a.GetName()).IsValid():
2188 aa_local_lddt[cname][rnum][a.GetName()] =
None
2189 elif trg_r
is not None and trg_r.FindAtom(a.GetName()).IsValid()
and \
2190 mdl_r
is not None and not mdl_r.FindAtom(a.GetName()).IsValid():
2192 aa_local_lddt[cname][rnum][a.GetName()] = 0.0
2195 aa_local_lddt[cname][rnum][a.GetName()] =
None
2199 if mdl_res.IsValid():
2202 local_lddt[cname][rnum] =
None
2204 aa_local_lddt[cname][rnum][a.GetName()] =
None
2211 if cname
in flat_mapping:
2212 for col
in alns[cname]:
2213 if col[0] !=
'-' and col[1] !=
'-':
2214 if col.GetResidue(1).number == r.number:
2215 trg_r = col.GetResidue(0)
2217 if trg_r
is not None:
2218 trg_cname = trg_r.GetChain().GetName()
2219 trg_rnum = trg_r.GetNumber()
2226 local_lddt[cname][rnum] =
None
2228 aa_local_lddt[cname][rnum][a.GetName()] =
None
2230 local_lddt[cname][rnum] = 0.0
2232 if trg_r.FindAtom(a.GetName()).IsValid():
2233 aa_local_lddt[cname][rnum][a.GetName()] = 0.0
2235 aa_local_lddt[cname][rnum][a.GetName()] =
None
2237 self.
_lddt_lddt = lddt_score
2241 def _compute_bb_lddt(self):
2242 LogScript(
"Computing backbone lDDT")
2245 for aln
in self.
alnaln:
2246 mdl_seq = aln.GetSequence(1)
2247 alns[mdl_seq.name] = aln
2250 flat_mapping = self.
mappingmapping.GetFlatMapping(mdl_as_key=
True)
2251 lddt_chain_mapping = dict()
2252 for mdl_ch, trg_ch
in flat_mapping.items():
2254 lddt_chain_mapping[mdl_ch] = trg_ch
2257 chain_mapping = lddt_chain_mapping,
2258 residue_mapping = alns,
2259 check_resnames=
False,
2260 local_lddt_prop=
"bb_lddt",
2263 for r
in self.
modelmodel.residues:
2264 cname = r.GetChain().GetName()
2265 if cname
not in local_lddt:
2266 local_lddt[cname] = dict()
2267 if r.HasProp(
"bb_lddt"):
2268 score = round(r.GetFloatProp(
"bb_lddt"), 3)
2269 local_lddt[cname][r.GetNumber()] = score
2273 local_lddt[cname][r.GetNumber()] =
None
2278 def _compute_ilddt(self):
2279 LogScript(
"Computing all-atom ilDDT")
2281 flat_mapping = self.
mappingmapping.GetFlatMapping(mdl_as_key=
True)
2285 for aln
in self.
alnaln:
2286 mdl_seq = aln.GetSequence(1)
2287 alns[mdl_seq.name] = aln
2288 lddt_chain_mapping = dict()
2289 for mdl_ch, trg_ch
in flat_mapping.items():
2291 lddt_chain_mapping[mdl_ch] = trg_ch
2293 chain_mapping = lddt_chain_mapping,
2294 residue_mapping = alns,
2295 check_resnames=
False,
2296 local_lddt_prop=
"lddt",
2298 no_intrachain=
True)[0]
2302 mdl_seq = aln.GetSequence(1)
2303 alns[mdl_seq.name] = aln
2304 lddt_chain_mapping = dict()
2305 for mdl_ch, trg_ch
in flat_mapping.items():
2307 lddt_chain_mapping[mdl_ch] = trg_ch
2309 chain_mapping = lddt_chain_mapping,
2310 residue_mapping = alns,
2311 check_resnames=
False,
2312 local_lddt_prop=
"lddt",
2314 no_intrachain=
True)[0]
2317 def _compute_qs(self):
2318 LogScript(
"Computing global QS-score")
2319 qs_score_result = self.
qs_scorerqs_scorer.Score(self.
mappingmapping.mapping)
2320 self.
_qs_global_qs_global = qs_score_result.QS_global
2321 self.
_qs_best_qs_best = qs_score_result.QS_best
2323 def _compute_per_interface_qs_scores(self):
2324 LogScript(
"Computing per-interface QS-score")
2329 trg_ch1 = interface[0]
2330 trg_ch2 = interface[1]
2331 mdl_ch1 = interface[2]
2332 mdl_ch2 = interface[3]
2333 qs_res = self.
qs_scorerqs_scorer.ScoreInterface(trg_ch1, trg_ch2,
2338 def _compute_ics_scores(self):
2339 LogScript(
"Computing ICS scores")
2342 self.
_ics_recall_ics_recall = contact_scorer_res.recall
2343 self.
_ics_ics = contact_scorer_res.ics
2347 flat_mapping = self.
mappingmapping.GetFlatMapping()
2349 trg_ch1 = trg_int[0]
2350 trg_ch2 = trg_int[1]
2351 if trg_ch1
in flat_mapping
and trg_ch2
in flat_mapping:
2352 mdl_ch1 = flat_mapping[trg_ch1]
2353 mdl_ch2 = flat_mapping[trg_ch2]
2354 res = self.
contact_scorercontact_scorer.ScoreICSInterface(trg_ch1, trg_ch2,
2364 def _trim_model(self):
2365 trimmed_mdl = mol.CreateEntityFromView(self.
mappingmapping.model,
True)
2366 trimmed_aln = dict()
2368 for trg_cname, mdl_cname
in self.
mappingmapping.GetFlatMapping().items():
2369 mdl_ch = trimmed_mdl.FindChain(mdl_cname)
2370 aln = self.
mappingmapping.alns[(trg_cname, mdl_cname)]
2373 assert(mdl_ch.GetResidueCount() == \
2374 len(aln.GetSequence(1).GetGaplessString()))
2376 mdl_residues = mdl_ch.residues
2378 aligned_mdl_seq = [
'-'] * aln.GetLength()
2379 for col_idx, col
in enumerate(aln):
2380 if col[0] !=
'-' and col[1] !=
'-':
2381 mdl_res = mdl_residues[mdl_res_idx]
2382 mdl_res.SetIntProp(
"aligned", 1)
2383 aligned_mdl_seq[col_idx] = col[1]
2386 aligned_mdl_seq =
''.join(aligned_mdl_seq)
2387 aligned_mdl_seq = seq.CreateSequence(mdl_cname, aligned_mdl_seq)
2389 new_aln = seq.CreateAlignment()
2390 new_aln.AddSequence(aln.GetSequence(0))
2391 new_aln.AddSequence(aligned_mdl_seq)
2392 trimmed_aln[(trg_cname, mdl_cname)] = new_aln
2394 self.
_trimmed_model_trimmed_model = trimmed_mdl.Select(
"graligned:0=1")
2397 def _compute_ics_scores_trimmed(self):
2398 LogScript(
"Computing ICS scores trimmed")
2412 flat_mapping = self.
mappingmapping.GetFlatMapping()
2414 trg_ch1 = trg_int[0]
2415 trg_ch2 = trg_int[1]
2416 if trg_ch1
in flat_mapping
and trg_ch2
in flat_mapping:
2417 mdl_ch1 = flat_mapping[trg_ch1]
2418 mdl_ch2 = flat_mapping[trg_ch2]
2429 def _compute_ips_scores(self):
2430 LogScript(
"Computing IPS scores")
2433 self.
_ips_recall_ips_recall = contact_scorer_res.recall
2434 self.
_ips_ips = contact_scorer_res.ips
2439 flat_mapping = self.
mappingmapping.GetFlatMapping()
2441 trg_ch1 = trg_int[0]
2442 trg_ch2 = trg_int[1]
2443 if trg_ch1
in flat_mapping
and trg_ch2
in flat_mapping:
2444 mdl_ch1 = flat_mapping[trg_ch1]
2445 mdl_ch2 = flat_mapping[trg_ch2]
2446 res = self.
contact_scorercontact_scorer.ScoreIPSInterface(trg_ch1, trg_ch2,
2456 def _compute_ips_scores_trimmed(self):
2457 LogScript(
"Computing IPS scores trimmed")
2470 flat_mapping = self.
mappingmapping.GetFlatMapping()
2472 trg_ch1 = trg_int[0]
2473 trg_ch2 = trg_int[1]
2474 if trg_ch1
in flat_mapping
and trg_ch2
in flat_mapping:
2475 mdl_ch1 = flat_mapping[trg_ch1]
2476 mdl_ch2 = flat_mapping[trg_ch2]
2487 def _compute_dockq_scores(self):
2488 LogScript(
"Computing DockQ")
2491 raise RuntimeError(
"Cannot compute DockQ for reference structures "
2492 "with nucleotide chains if dockq_capri_peptide "
2497 self.
_fnat_fnat = list()
2499 self.
_irmsd_irmsd = list()
2500 self.
_lrmsd_lrmsd = list()
2501 self.
_nnat_nnat = list()
2502 self.
_nmdl_nmdl = list()
2505 for aln
in self.
alnaln:
2506 trg_s = aln.GetSequence(0)
2507 mdl_s = aln.GetSequence(1)
2508 dockq_alns[(trg_s.GetName(), mdl_s.GetName())] = aln
2511 trg_ch1 = interface[0]
2512 trg_ch2 = interface[1]
2513 mdl_ch1 = interface[2]
2514 mdl_ch2 = interface[3]
2515 aln1 = dockq_alns[(trg_ch1, mdl_ch1)]
2516 aln2 = dockq_alns[(trg_ch2, mdl_ch2)]
2518 res = dockq.DockQ(self.
modelmodel, self.
targettarget, mdl_ch1, mdl_ch2,
2519 trg_ch1, trg_ch2, ch1_aln=aln1,
2520 ch2_aln=aln2, contact_dist_thresh = 4.0,
2521 interface_dist_thresh=8.0,
2522 interface_cb =
True)
2524 res = dockq.DockQ(self.
modelmodel, self.
targettarget, mdl_ch1, mdl_ch2,
2525 trg_ch1, trg_ch2, ch1_aln=aln1,
2528 self.
_fnat_fnat.append(res[
"fnat"])
2529 self.
_fnonnat_fnonnat.append(res[
"fnonnat"])
2530 self.
_irmsd_irmsd.append(res[
"irmsd"])
2531 self.
_lrmsd_lrmsd.append(res[
"lrmsd"])
2533 self.
_nnat_nnat.append(res[
"nnat"])
2534 self.
_nmdl_nmdl.append(res[
"nmdl"])
2539 not_covered_counts = list()
2540 proc_trg_interfaces = set([(x[0], x[1])
for x
in self.
dockq_interfacesdockq_interfaces])
2542 if interface
not in proc_trg_interfaces:
2546 trg_ch1 = interface[0]
2547 trg_ch2 = interface[1]
2550 res = dockq.DockQ(self.
targettarget, self.
targettarget,
2551 trg_ch1, trg_ch2, trg_ch1, trg_ch2,
2552 contact_dist_thresh = 4.0,
2553 interface_dist_thresh=8.0,
2554 interface_cb =
True)
2556 res = dockq.DockQ(self.
targettarget, self.
targettarget,
2557 trg_ch1, trg_ch2, trg_ch1, trg_ch2)
2559 not_covered_counts.append(res[
"nnat"])
2566 weights = np.array(self.
_nnat_nnat)
2571 self.
_dockq_wave_dockq_wave = np.sum(np.multiply(weights/np.sum(weights), scores))
2572 scores = np.append(scores, [0.0]*len(not_covered_counts))
2573 weights = np.append(weights, not_covered_counts)
2578 self.
_dockq_wave_full_dockq_wave_full = np.sum(np.multiply(weights/np.sum(weights),
2581 def _extract_mapped_pos(self):
2585 processed_trg_chains = set()
2586 for trg_ch, mdl_ch
in self.
mappingmapping.GetFlatMapping().items():
2587 processed_trg_chains.add(trg_ch)
2588 aln = self.
mappingmapping.alns[(trg_ch, mdl_ch)]
2590 if col[0] !=
'-' and col[1] !=
'-':
2591 trg_res = col.GetResidue(0)
2592 mdl_res = col.GetResidue(1)
2593 trg_at = trg_res.FindAtom(
"CA")
2594 mdl_at = mdl_res.FindAtom(
"CA")
2595 if not trg_at.IsValid():
2596 trg_at = trg_res.FindAtom(
"C3'")
2597 if not mdl_at.IsValid():
2598 mdl_at = mdl_res.FindAtom(
"C3'")
2604 for ch
in self.
mappingmapping.target.chains:
2605 if ch.GetName()
not in processed_trg_chains:
2608 def _extract_mapped_pos_full_bb(self):
2611 exp_pep_atoms = [
"N",
"CA",
"C"]
2612 exp_nuc_atoms = [
"\"O5'\"",
"\"C5'\"",
"\"C4'\"",
"\"C3'\"",
"\"O3'\""]
2613 trg_pep_chains = [s.GetName()
for s
in self.
chain_mapperchain_mapper.polypep_seqs]
2614 trg_nuc_chains = [s.GetName()
for s
in self.
chain_mapperchain_mapper.polynuc_seqs]
2615 for trg_ch, mdl_ch
in self.
mappingmapping.GetFlatMapping().items():
2616 aln = self.
mappingmapping.alns[(trg_ch, mdl_ch)]
2617 trg_ch = aln.GetSequence(0).GetName()
2618 if trg_ch
in trg_pep_chains:
2619 exp_atoms = exp_pep_atoms
2620 elif trg_ch
in trg_nuc_chains:
2621 exp_atoms = exp_nuc_atoms
2624 raise RuntimeError(
"Unexpected error - contact OST developer")
2626 if col[0] !=
'-' and col[1] !=
'-':
2627 trg_res = col.GetResidue(0)
2628 mdl_res = col.GetResidue(1)
2629 for aname
in exp_atoms:
2630 trg_at = trg_res.FindAtom(aname)
2631 mdl_at = mdl_res.FindAtom(aname)
2632 if not (trg_at.IsValid()
and mdl_at.IsValid()):
2634 raise RuntimeError(
"Unexpected error - contact OST "
2640 def _extract_rigid_mapped_pos(self):
2644 processed_trg_chains = set()
2645 for trg_ch, mdl_ch
in self.
rigid_mappingrigid_mapping.GetFlatMapping().items():
2646 processed_trg_chains.add(trg_ch)
2647 aln = self.
rigid_mappingrigid_mapping.alns[(trg_ch, mdl_ch)]
2649 if col[0] !=
'-' and col[1] !=
'-':
2650 trg_res = col.GetResidue(0)
2651 mdl_res = col.GetResidue(1)
2652 trg_at = trg_res.FindAtom(
"CA")
2653 mdl_at = mdl_res.FindAtom(
"CA")
2654 if not trg_at.IsValid():
2655 trg_at = trg_res.FindAtom(
"C3'")
2656 if not mdl_at.IsValid():
2657 mdl_at = mdl_res.FindAtom(
"C3'")
2664 if ch.GetName()
not in processed_trg_chains:
2667 def _extract_rigid_mapped_pos_full_bb(self):
2670 exp_pep_atoms = [
"N",
"CA",
"C"]
2671 exp_nuc_atoms = [
"\"O5'\"",
"\"C5'\"",
"\"C4'\"",
"\"C3'\"",
"\"O3'\""]
2672 trg_pep_chains = [s.GetName()
for s
in self.
chain_mapperchain_mapper.polypep_seqs]
2673 trg_nuc_chains = [s.GetName()
for s
in self.
chain_mapperchain_mapper.polynuc_seqs]
2674 for trg_ch, mdl_ch
in self.
rigid_mappingrigid_mapping.GetFlatMapping().items():
2675 aln = self.
mappingmapping.alns[(trg_ch, mdl_ch)]
2676 trg_ch = aln.GetSequence(0).GetName()
2677 if trg_ch
in trg_pep_chains:
2678 exp_atoms = exp_pep_atoms
2679 elif trg_ch
in trg_nuc_chains:
2680 exp_atoms = exp_nuc_atoms
2683 raise RuntimeError(
"Unexpected error - contact OST developer")
2685 if col[0] !=
'-' and col[1] !=
'-':
2686 trg_res = col.GetResidue(0)
2687 mdl_res = col.GetResidue(1)
2688 for aname
in exp_atoms:
2689 trg_at = trg_res.FindAtom(aname)
2690 mdl_at = mdl_res.FindAtom(aname)
2691 if not (trg_at.IsValid()
and mdl_at.IsValid()):
2693 raise RuntimeError(
"Unexpected error - contact OST "
2698 def _compute_cad_score(self):
2700 raise RuntimeError(
"CAD score computations rely on residue numbers "
2701 "that are consistent between target and model "
2702 "chains, i.e. only work if resnum_alignments "
2703 "is True at Scorer construction.")
2705 LogScript(
"Computing CAD score")
2707 settings.Locate(
"voronota-cadscore",
2709 except Exception
as e:
2710 raise RuntimeError(
"voronota-cadscore must be in PATH for CAD "
2711 "score scoring")
from e
2712 cad_bin_dir = os.path.dirname(cad_score_exec)
2713 m = self.
mappingmapping.GetFlatMapping(mdl_as_key=
True)
2714 cad_result = cadscore.CADScore(self.
modelmodel, self.
targettarget,
2718 cad_bin_path=cad_bin_dir,
2722 for r
in self.
modelmodel.residues:
2723 cname = r.GetChain().GetName()
2724 if cname
not in local_cad:
2725 local_cad[cname] = dict()
2726 if r.HasProp(
"localcad"):
2727 score = round(r.GetFloatProp(
"localcad"), 3)
2728 local_cad[cname][r.GetNumber()] = score
2730 local_cad[cname][r.GetNumber()] =
None
2732 self.
_cad_score_cad_score = cad_result.globalAA
2735 def _get_repr_view(self, ent):
2736 """ Returns view with representative peptide atoms => CB, CA for GLY
2738 Ensures that each residue has exactly one atom with assertions
2740 :param ent: Entity for which you want the representative view
2741 :param ent: :class:`ost.mol.EntityHandle`/:class:`ost.mol.EntityView`
2742 :returns: :class:`ost.mol.EntityView` derived from ent
2744 repr_ent = ent.Select(
"(aname=\"CB\" or (rname=\"GLY\" and aname=\"CA\"))")
2745 for r
in repr_ent.residues:
2746 assert(len(r.atoms) == 1)
2749 def _get_interface_residues(self, ent):
2750 """ Get interface residues
2752 Thats all residues having a contact with at least one residue from another
2753 chain (CB-CB distance <= 8A, CA in case of Glycine)
2755 :param ent: Model for which to extract interface residues
2756 :type ent: :class:`ost.mol.EntityView`
2757 :returns: :class:`dict` with chain names as key and and :class:`list`
2758 with residue numbers of the respective interface residues.
2762 result = {ch.GetName(): list()
for ch
in ent.chains}
2763 for ch
in ent.chains:
2764 cname = ch.GetName()
2765 sel = repr_ent.Select(f
"(cname={mol.QueryQuoteName(cname)} and 8 <> [cname!={mol.QueryQuoteName(cname)}])")
2766 result[cname] = [r.GetNumber()
for r
in sel.residues]
2769 def _do_stereochecks(self):
2770 """ Perform stereochemistry checks on model and target
2772 LogInfo(
"Performing stereochemistry checks on model and target")
2773 data = stereochemistry.GetDefaultStereoData()
2774 l_data = stereochemistry.GetDefaultStereoLinkData()
2776 a, b, c, d = stereochemistry.StereoCheck(self.
modelmodel, stereo_data = data,
2777 stereo_link_data = l_data)
2783 a, b, c, d = stereochemistry.StereoCheck(self.
targettarget, stereo_data = data,
2784 stereo_link_data = l_data)
2791 def _get_interface_patches(self, mdl_ch, mdl_rnum):
2792 """ Select interface patches representative for specified residue
2794 The patches for specified residue r in chain c are selected as follows:
2796 * mdl_patch_one: All residues in c with CB (CA for GLY) positions within 8A
2797 of r and within 12A of residues from any other chain.
2798 * mdl_patch_two: Closest residue x to r in any other chain gets identified.
2799 Patch is then constructed by selecting all residues from any other chain
2800 within 8A of x and within 12A from any residue in c.
2801 * trg_patch_one: Chain name and residue number based mapping from
2803 * trg_patch_two: Chain name and residue number based mapping from
2806 :param mdl_ch: Name of chain in *self.model* of residue of interest
2807 :type mdl_ch: :class:`str`
2808 :param mdl_rnum: Residue number of residue of interest
2809 :type mdl_rnum: :class:`ost.mol.ResNum`
2810 :returns: Tuple with 5 elements: 1) :class:`bool` flag whether all residues
2811 in *mdl* patches are covered in *trg* 2) mtl_patch_one
2812 3) mdl_patch_two 4) trg_patch_one 5) trg_patch_two
2818 r = self.
modelmodel.FindResidue(mdl_ch, mdl_rnum)
2820 raise RuntimeError(f
"Cannot find residue {mdl_rnum} in chain {mdl_ch}")
2821 if r.GetName() ==
"GLY":
2822 at = r.FindAtom(
"CA")
2824 at = r.FindAtom(
"CB")
2825 if not at.IsValid():
2826 raise RuntimeError(
"Cannot find interface views for res without CB/CA")
2833 q1 = f
"(cname={mol.QueryQuoteName(mdl_ch)} and 8 <> {{{r_pos[0]},{r_pos[1]},{r_pos[2]}}})"
2835 q2 = f
"(12 <> [cname!={mol.QueryQuoteName(mdl_ch)}])"
2836 mdl_patch_one = self.
modelmodel.CreateEmptyView()
2837 sel = repr_mdl.Select(
" and ".join([q1, q2]))
2838 for r
in sel.residues:
2839 mdl_r = self.
modelmodel.FindResidue(r.GetChain().GetName(), r.GetNumber())
2840 mdl_patch_one.AddResidue(mdl_r, mol.ViewAddFlag.INCLUDE_ALL)
2846 sel = repr_mdl.Select(f
"(cname!={mol.QueryQuoteName(mdl_ch)})")
2847 close_stuff = sel.FindWithin(r_pos, 8)
2850 for close_at
in close_stuff:
2853 min_pos = close_at.GetPos()
2857 q1 = f
"(cname!={mol.QueryQuoteName(mdl_ch)} and 8 <> {{{min_pos[0]},{min_pos[1]},{min_pos[2]}}})"
2859 q2 = f
"(12 <> [cname={mol.QueryQuoteName(mdl_ch)}])"
2860 mdl_patch_two = self.
modelmodel.CreateEmptyView()
2861 sel = repr_mdl.Select(
" and ".join([q1, q2]))
2862 for r
in sel.residues:
2863 mdl_r = self.
modelmodel.FindResidue(r.GetChain().GetName(), r.GetNumber())
2864 mdl_patch_two.AddResidue(mdl_r, mol.ViewAddFlag.INCLUDE_ALL)
2867 flat_mapping = self.
mappingmapping.GetFlatMapping(mdl_as_key=
True)
2868 full_trg_coverage =
True
2869 trg_patch_one = self.
targettarget.CreateEmptyView()
2870 for r
in mdl_patch_one.residues:
2872 mdl_cname = r.GetChain().GetName()
2873 if mdl_cname
in flat_mapping:
2874 aln = self.
mappingmapping.alns[(flat_mapping[mdl_cname], mdl_cname)]
2876 if col[0] !=
'-' and col[1] !=
'-':
2877 if col.GetResidue(1).GetNumber() == r.GetNumber():
2878 trg_r = col.GetResidue(0)
2880 if trg_r
is not None:
2881 trg_patch_one.AddResidue(trg_r.handle,
2882 mol.ViewAddFlag.INCLUDE_ALL)
2884 full_trg_coverage =
False
2886 trg_patch_two = self.
targettarget.CreateEmptyView()
2887 for r
in mdl_patch_two.residues:
2889 mdl_cname = r.GetChain().GetName()
2890 if mdl_cname
in flat_mapping:
2891 aln = self.
mappingmapping.alns[(flat_mapping[mdl_cname], mdl_cname)]
2893 if col[0] !=
'-' and col[1] !=
'-':
2894 if col.GetResidue(1).GetNumber() == r.GetNumber():
2895 trg_r = col.GetResidue(0)
2897 if trg_r
is not None:
2898 trg_patch_two.AddResidue(trg_r.handle,
2899 mol.ViewAddFlag.INCLUDE_ALL)
2901 full_trg_coverage =
False
2903 return (full_trg_coverage, mdl_patch_one, mdl_patch_two,
2904 trg_patch_one, trg_patch_two)
2906 def _compute_patchqs_scores(self):
2907 LogScript(
"Computing patch QS-scores")
2913 r = self.
modelmodel.FindResidue(cname, rnum)
2914 if r.IsValid()
and r.GetChemType() == mol.ChemType.AMINOACIDS:
2915 full_trg_coverage, mdl_patch_one, mdl_patch_two, \
2916 trg_patch_one, trg_patch_two = \
2918 if full_trg_coverage:
2919 score = self.
_patchqs_patchqs(mdl_patch_one, mdl_patch_two,
2920 trg_patch_one, trg_patch_two)
2921 scores.append(score)
2924 def _compute_patchdockq_scores(self):
2925 LogScript(
"Computing patch DockQ scores")
2931 r = self.
modelmodel.FindResidue(cname, rnum)
2932 if r.IsValid()
and r.GetChemType() == mol.ChemType.AMINOACIDS:
2933 full_trg_coverage, mdl_patch_one, mdl_patch_two, \
2934 trg_patch_one, trg_patch_two = \
2936 if full_trg_coverage:
2937 score = self.
_patchdockq_patchdockq(mdl_patch_one, mdl_patch_two,
2938 trg_patch_one, trg_patch_two)
2939 scores.append(score)
2942 def _patchqs(self, mdl_patch_one, mdl_patch_two, trg_patch_one, trg_patch_two):
2943 """ Score interface residue patches with QS-score
2945 In detail: Construct two entities with two chains each. First chain
2946 consists of residues from <x>_patch_one and second chain consists of
2947 <x>_patch_two. The returned score is the QS-score between the two
2950 :param mdl_patch_one: Interface patch representing scored residue
2951 :type mdl_patch_one: :class:`ost.mol.EntityView`
2952 :param mdl_patch_two: Interface patch representing scored residue
2953 :type mdl_patch_two: :class:`ost.mol.EntityView`
2954 :param trg_patch_one: Interface patch representing scored residue
2955 :type trg_patch_one: :class:`ost.mol.EntityView`
2956 :param trg_patch_two: Interface patch representing scored residue
2957 :type trg_patch_two: :class:`ost.mol.EntityView`
2958 :returns: PatchQS score
2963 alnA = seq.CreateAlignment()
2964 s =
''.join([r.one_letter_code
for r
in mdl_patch_one.residues])
2965 alnA.AddSequence(seq.CreateSequence(
"A", s))
2966 s =
''.join([r.one_letter_code
for r
in trg_patch_one.residues])
2967 alnA.AddSequence(seq.CreateSequence(
"A", s))
2969 alnB = seq.CreateAlignment()
2970 s =
''.join([r.one_letter_code
for r
in mdl_patch_two.residues])
2971 alnB.AddSequence(seq.CreateSequence(
"B", s))
2972 s =
''.join([r.one_letter_code
for r
in trg_patch_two.residues])
2973 alnB.AddSequence(seq.CreateSequence(
"B", s))
2974 alns = {(
"A",
"A"): alnA, (
"B",
"B"): alnB}
2976 scorer =
QSScorer(qs_ent_mdl, [[
"A"], [
"B"]], qs_ent_trg, alns)
2977 score_result = scorer.Score([[
"A"], [
"B"]])
2979 return score_result.QS_global
2981 def _patchdockq(self, mdl_patch_one, mdl_patch_two, trg_patch_one,
2983 """ Score interface residue patches with DockQ
2985 In detail: Construct two entities with two chains each. First chain
2986 consists of residues from <x>_patch_one and second chain consists of
2987 <x>_patch_two. The returned score is the QS-score between the two
2990 :param mdl_patch_one: Interface patch representing scored residue
2991 :type mdl_patch_one: :class:`ost.mol.EntityView`
2992 :param mdl_patch_two: Interface patch representing scored residue
2993 :type mdl_patch_two: :class:`ost.mol.EntityView`
2994 :param trg_patch_one: Interface patch representing scored residue
2995 :type trg_patch_one: :class:`ost.mol.EntityView`
2996 :param trg_patch_two: Interface patch representing scored residue
2997 :type trg_patch_two: :class:`ost.mol.EntityView`
2998 :returns: DockQ score
3002 dockq_result = dockq.DockQ(t, m,
"A",
"B",
"A",
"B")
3003 if dockq_result[
"nnat"] > 0:
3004 return dockq_result[
"DockQ"]
3007 def _qs_ent_from_patches(self, patch_one, patch_two):
3008 """ Constructs Entity with two chains named "A" and "B""
3010 Blindly adds all residues from *patch_one* to chain A and residues from
3011 patch_two to chain B.
3013 ent = mol.CreateEntity()
3015 added_ch = ed.InsertChain(
"A")
3016 for r
in patch_one.residues:
3017 added_r = ed.AppendResidue(added_ch, r.GetName())
3018 added_r.SetChemClass(str(r.GetChemClass()))
3020 ed.InsertAtom(added_r, a.handle)
3021 added_ch = ed.InsertChain(
"B")
3022 for r
in patch_two.residues:
3023 added_r = ed.AppendResidue(added_ch, r.GetName())
3024 added_r.SetChemClass(str(r.GetChemClass()))
3026 ed.InsertAtom(added_r, a.handle)
3029 def _construct_custom_mapping(self, mapping):
3030 """ constructs a full blown MappingResult object from simple dict
3032 :param mapping: mapping with trg chains as key and mdl ch as values
3033 :type mapping: :class:`dict`
3035 LogInfo(
"Setting custom chain mapping")
3038 chem_mapping, chem_group_alns, mdl = \
3039 chain_mapper.GetChemMapping(self.
modelmodel)
3044 if len(mapping) != len(set(mapping.keys())):
3045 raise RuntimeError(f
"Expect unique trg chain names in mapping. Got "
3046 f
"{mapping.keys()}")
3047 if len(mapping) != len(set(mapping.values())):
3048 raise RuntimeError(f
"Expect unique mdl chain names in mapping. Got "
3049 f
"{mapping.values()}")
3051 trg_chains = set([ch.GetName()
for ch
in chain_mapper.target.chains])
3052 mdl_chains = set([ch.GetName()
for ch
in mdl.chains])
3053 for k,v
in mapping.items():
3054 if k
not in trg_chains:
3055 raise RuntimeError(f
"Target chain \"{k}\" is not available "
3056 f
"in target processed for chain mapping "
3058 if v
not in mdl_chains:
3059 raise RuntimeError(f
"Model chain \"{v}\" is not available "
3060 f
"in model processed for chain mapping "
3063 for trg_ch, mdl_ch
in mapping.items():
3064 trg_group_idx =
None
3065 mdl_group_idx =
None
3066 for idx, group
in enumerate(chain_mapper.chem_groups):
3070 for idx, group
in enumerate(chem_mapping):
3074 if trg_group_idx
is None or mdl_group_idx
is None:
3075 raise RuntimeError(
"Could not establish a valid chem grouping "
3076 "of chain names provided in custom mapping.")
3078 if trg_group_idx != mdl_group_idx:
3079 raise RuntimeError(f
"Chem group mismatch in custom mapping: "
3080 f
"target chain \"{trg_ch}\" groups with the "
3081 f
"following chemically equivalent target "
3083 f
"{chain_mapper.chem_groups[trg_group_idx]} "
3084 f
"but model chain \"{mdl_ch}\" maps to the "
3085 f
"following target chains: "
3086 f
"{chain_mapper.chem_groups[mdl_group_idx]}")
3088 pairs = set([(trg_ch, mdl_ch)
for trg_ch, mdl_ch
in mapping.items()])
3090 chain_mapping._GetRefMdlAlns(chain_mapper.chem_groups,
3091 chain_mapper.chem_group_alignments,
3097 final_mapping = list()
3098 for ref_chains
in chain_mapper.chem_groups:
3099 mapped_mdl_chains = list()
3100 for ref_ch
in ref_chains:
3101 if ref_ch
in mapping:
3102 mapped_mdl_chains.append(mapping[ref_ch])
3104 mapped_mdl_chains.append(
None)
3105 final_mapping.append(mapped_mdl_chains)
3108 for ref_group, mdl_group
in zip(chain_mapper.chem_groups,
3110 for ref_ch, mdl_ch
in zip(ref_group, mdl_group):
3111 if ref_ch
is not None and mdl_ch
is not None:
3112 aln = ref_mdl_alns[(ref_ch, mdl_ch)]
3113 trg_view = chain_mapper.target.Select(f
"cname={mol.QueryQuoteName(ref_ch)}")
3114 mdl_view = mdl.Select(f
"cname={mol.QueryQuoteName(mdl_ch)}")
3115 aln.AttachView(0, trg_view)
3116 aln.AttachView(1, mdl_view)
3117 alns[(ref_ch, mdl_ch)] = aln
3120 chain_mapper.chem_groups,
3122 final_mapping, alns)
3124 def _compute_tmscore(self):
3127 LogScript(
"Computing TM-score with built-in USalign")
3129 flat_mapping = self.
rigid_mappingrigid_mapping.GetFlatMapping()
3130 LogInfo(
"Overriding TM-score chain mapping")
3131 res = res = bindings.WrappedMMAlign(self.
modelmodel, self.
targettarget,
3132 mapping=flat_mapping)
3134 res = bindings.WrappedMMAlign(self.
modelmodel, self.
targettarget)
3136 LogScript(
"Computing TM-score with USalign exectuable")
3138 LogInfo(
"Overriding TM-score chain mapping")
3139 flat_mapping = self.
rigid_mappingrigid_mapping.GetFlatMapping()
3140 res = tmtools.USAlign(self.
modelmodel, self.
targettarget,
3142 custom_chain_mapping = flat_mapping)
3144 res = tmtools.USAlign(self.
modelmodel, self.
targettarget,
3149 for a,b
in zip(res.ent1_mapped_chains, res.ent2_mapped_chains):
3152 def _resnum_connect(self, ent):
3154 for ch
in ent.chains:
3155 res_list = ch.residues
3156 for i
in range(len(res_list) - 1):
3159 if ra.GetNumber().GetNum() + 1 == rb.GetNumber().GetNum():
3160 if ra.IsPeptideLinking()
and rb.IsPeptideLinking():
3161 c = ra.FindAtom(
"C")
3162 n = rb.FindAtom(
"N")
3163 if c.IsValid()
and n.IsValid()
and not mol.BondExists(c, n):
3165 ed = ent.EditXCS(mol.BUFFERED_EDIT)
3167 elif ra.IsNucleotideLinking()
and rb.IsNucleotideLinking():
3168 o = ra.FindAtom(
"O3'")
3169 p = rb.FindAtom(
"P")
3170 if o.IsValid()
and p.IsValid()
and not mol.BondExists(o, p):
3172 ed = ent.EditXCS(mol.BUFFERED_EDIT)
3177 __all__ = (
'lDDTBSScorer',
'Scorer',)
def target_interface_residues(self)
_per_interface_ips_trimmed
def _patchqs(self, mdl_patch_one, mdl_patch_two, trg_patch_one, trg_patch_two)
_per_interface_ips_recall
def per_interface_ics(self)
def target_bad_angles(self)
def _get_repr_view(self, ent)
def qs_target_interfaces(self)
def rigid_mapped_model_pos(self)
def _do_stereochecks(self)
def _compute_per_interface_qs_scores(self)
def model_bad_angles(self)
_rigid_mapped_target_pos_full_bb
_per_interface_ics_recall
_per_interface_ips_precision
def model_interface_residues(self)
_contact_target_interfaces
def stereochecked_target(self)
def per_interface_ips_precision(self)
_mapped_model_pos_full_bb
def dockq_wave_full(self)
def mapped_model_pos(self)
def rigid_mapped_target_pos_full_bb(self)
_per_interface_ics_precision_trimmed
def _compute_dockq_scores(self)
def _aln_helper(self, target, model)
def model_bad_bonds(self)
def per_interface_ics_recall(self)
def ics_precision_trimmed(self)
def _compute_bb_lddt(self)
def transformed_mapped_model_pos(self)
_per_interface_ips_recall_trimmed
_per_interface_ics_recall_trimmed
def usalign_mapping(self)
def native_contacts(self)
def local_cad_score(self)
def _compute_patchqs_scores(self)
def per_interface_ips(self)
def _compute_stereochecked_aln(self)
_rigid_transformed_mapped_model_pos
def rigid_transformed_mapped_model_pos(self)
def trimmed_model_contacts(self)
def _compute_ips_scores_trimmed(self)
def n_target_not_mapped(self)
def per_interface_ips_trimmed(self)
def contact_model_interfaces(self)
_transformed_mapped_model_pos
def _patchdockq(self, mdl_patch_one, mdl_patch_two, trg_patch_one, trg_patch_two)
_rigid_mapped_model_pos_full_bb
def dockq_interfaces(self)
def per_interface_ips_recall(self)
def per_interface_ics_precision(self)
def _compute_ics_scores(self)
def _resnum_connect(self, ent)
_rigid_n_target_not_mapped
def _get_interface_residues(self, ent)
def ips_precision_trimmed(self)
def rigid_transform(self)
_mapped_target_pos_full_bb
def ics_recall_trimmed(self)
def mapped_target_pos_full_bb(self)
def _compute_ips_scores(self)
def _extract_mapped_pos(self)
def __init__(self, model, target, resnum_alignments=False, molck_settings=None, cad_score_exec=None, custom_mapping=None, custom_rigid_mapping=None, usalign_exec=None, lddt_no_stereochecks=False, n_max_naive=40320, oum=False, min_pep_length=6, min_nuc_length=4, lddt_add_mdl_contacts=False, dockq_capri_peptide=False, lddt_symmetry_settings=None, lddt_inclusion_radius=15.0)
def rigid_mapped_target_pos(self)
def _compute_tmscore(self)
def dockq_target_interfaces(self)
def stereochecked_aln(self)
def trimmed_contact_scorer(self)
def per_interface_ics_recall_trimmed(self)
_per_interface_ics_precision
def _compute_patchdockq_scores(self)
def _get_interface_patches(self, mdl_ch, mdl_rnum)
def _extract_rigid_mapped_pos(self)
_contact_model_interfaces
_target_interface_residues
def contact_target_interfaces(self)
def per_interface_qs_global(self)
def _qs_ent_from_patches(self, patch_one, patch_two)
def rigid_mapped_model_pos_full_bb(self)
def stereochecked_model(self)
_model_interface_residues
def _compute_pepnuc_aln(self)
def per_interface_qs_best(self)
def _extract_rigid_mapped_pos_full_bb(self)
def per_interface_ips_precision_trimmed(self)
def _compute_ics_scores_trimmed(self)
_per_interface_ips_precision_trimmed
def rigid_n_target_not_mapped(self)
def per_interface_ics_precision_trimmed(self)
def per_interface_ips_recall_trimmed(self)
def mapped_target_pos(self)
def mapped_model_pos_full_bb(self)
def target_bad_bonds(self)
def qs_model_interfaces(self)
def _compute_cad_score(self)
def ips_recall_trimmed(self)
def per_interface_ics_trimmed(self)
_per_interface_ics_trimmed
def _construct_custom_mapping(self, mapping)
def _extract_mapped_pos_full_bb(self)
def ScoreBS(self, ligand, radius=4.0, lddt_radius=10.0)
def __init__(self, reference, model, residue_number_alignment=False)
Real DLLEXPORT_OST_GEOM Distance(const Line2 &l, const Vec2 &v)
void Molck(ost::mol::EntityHandle &ent, ost::conop::CompoundLibPtr lib, const MolckSettings &settings, bool prune=true)
void GDT(const geom::Vec3List &mdl_pos, const geom::Vec3List &ref_pos, int window_size, int max_windows, Real distance_thresh, int &n_superposed, geom::Mat4 &transform)