Skip to content

Commit

Permalink
Fix bug in TM-score computation for short chains
Browse files Browse the repository at this point in the history
  • Loading branch information
padix-key committed Jan 17, 2025
1 parent 8138e0e commit c0a1b93
Show file tree
Hide file tree
Showing 2 changed files with 24 additions and 3 deletions.
10 changes: 7 additions & 3 deletions src/biotite/structure/tm.py
Original file line number Diff line number Diff line change
Expand Up @@ -533,8 +533,7 @@ def _tm_score(distance, reference_length):
tm_score : float or ndarray
The TM score for the given distances.
"""
d0 = max(_d0(reference_length), _D0_MIN)
return 1 / (1 + (distance / d0) ** 2)
return 1 / (1 + (distance / _d0(reference_length)) ** 2)


def _d0(reference_length):
Expand All @@ -552,7 +551,12 @@ def _d0(reference_length):
The :math:`d_0` threshold.
"""
# Constants taken from Zhang2004
return 1.24 * (reference_length - 15) ** (1 / 3) - 1.8
return max(
# Avoid complex solutions -> clip to positive values
# For short sequence lengths _D0_MIN takes precedence anyway
1.24 * max((reference_length - 15), 0) ** (1 / 3) - 1.8,
_D0_MIN,
)


def _get_reference_length(user_parameter, reference_length, subject_length):
Expand Down
17 changes: 17 additions & 0 deletions tests/structure/test_tm.py
Original file line number Diff line number Diff line change
Expand Up @@ -158,6 +158,23 @@ def test_superimpose_consistency(fixed_pdb_id, mobile_pdb_id, ref_tm_score):
)


def test_short_structure():
"""
Check that the :func:`tm_score()` and :func:`superimpose_structural_homologs()` work
for structures shorter than 15 residues, because in that case we got
``TypeError: '>' not supported between instances of 'float' and 'complex'``
before.
"""
pdbx_file = pdbx.BinaryCIFFile.read(Path(data_dir("structure")) / "1l2y.bcif")
atoms = pdbx.get_structure(pdbx_file, model=1)
atoms = atoms[atoms.res_id < 15]

_, _, fixed_i, mobile_i = struc.superimpose_structural_homologs(atoms, atoms)
struc.tm_score(atoms, atoms, fixed_i, mobile_i)


def _transform_random_affine(coord, rng):
coord = struc.translate(coord, rng.uniform(low=0, high=10, size=3))
coord = struc.rotate(coord, rng.uniform(low=0, high=2 * np.pi, size=3))
Expand Down

0 comments on commit c0a1b93

Please sign in to comment.