Hunt-Szymanski algorithm in Python
Computes LCS in O((r + n) log n) time, where r is the number of matching index pairs, which is fast when matches are sparse.
- File
hunt_szymanski.py, 55 lines- Level
- I0 not checked yet
- Toolchain
- Python 3.14, pytest 9
- Test framework
- pytest
- Run the tests
python3 -m pytest
# SPDX-License-Identifier: MIT
# Part of https://algorithm.tamnd.com/hunt-szymanski/. The authors are listed in prov/EEQDBV.jsonl.
from bisect import bisect_left
from collections.abc import Hashable, Sequence
from typing import TypeVar
T = TypeVar("T", bound=Hashable)
def match_lists(x: Sequence[T], y: Sequence[T]) -> list[list[int]]:
"""Return, for each i, the 1-based positions j with y[j] == x[i], in decreasing order."""
where: dict[T, list[int]] = {}
for j in range(len(y), 0, -1): # @line L1
where.setdefault(y[j - 1], []).append(j)
return [where.get(c, []) for c in x]
def hunt_szymanski(x: Sequence[T], y: Sequence[T]) -> list[T]:
"""Return one longest common subsequence of x and y."""
lists = match_lists(x, y)
thresh: list[int] = [0] # @line L2
link: list[tuple[int, int, object] | None] = [None]
for i in range(1, len(x) + 1): # @line L3
for j in lists[i - 1]: # @line L4
k = bisect_left(thresh, j) # @line L5
if k == len(thresh): # @line L6
thresh.append(j) # @line L7
link.append((i, j, link[k - 1]))
elif j < thresh[k]: # @line L6
thresh[k] = j # @line L7
link[k] = (i, j, link[k - 1])
z: list[T] = []
node = link[-1] # @line L8
while node is not None: # @line L9
i, j, node = node # type: ignore[misc]
z.append(x[i - 1])
z.reverse()
return z # @line L10
def hunt_szymanski_trace(x: Sequence[T], y: Sequence[T]) -> list[list[int]]:
"""Return the threshold row T[1 ..] after each item of x (for the example on the page)."""
lists = match_lists(x, y)
thresh = [0]
rows = []
for i in range(len(x)):
for j in lists[i]:
k = bisect_left(thresh, j)
if k == len(thresh):
thresh.append(j)
elif j < thresh[k]:
thresh[k] = j
rows.append(thresh[1:])
return rows
test_hunt_szymanski.py 56 lines
# SPDX-License-Identifier: MIT
# Part of https://algorithm.tamnd.com/hunt-szymanski/. The authors are listed in prov/EEQDBV.jsonl.
import random
from hunt_szymanski import hunt_szymanski, hunt_szymanski_trace, match_lists
def is_subsequence(z: list, x: list) -> bool:
it = iter(x)
return all(any(a == b for b in it) for a in z)
def table_lcs_length(x: str, y: str) -> int:
c = [[0] * (len(y) + 1) for _ in range(len(x) + 1)]
for i in range(1, len(x) + 1):
for j in range(1, len(y) + 1):
if x[i - 1] == y[j - 1]:
c[i][j] = c[i - 1][j - 1] + 1
else:
c[i][j] = max(c[i - 1][j], c[i][j - 1])
return c[len(x)][len(y)]
def test_example() -> None:
x, y = "ABCBDAB", "BDCABA"
assert match_lists(x, y)[0] == [6, 4]
assert sum(len(lst) for lst in match_lists(x, y)) == 12
rows = hunt_szymanski_trace(x, y)
assert rows == [[4], [1, 5], [1, 3], [1, 3, 5], [1, 2, 5], [1, 2, 4, 6], [1, 2, 4, 5]]
assert "".join(hunt_szymanski(x, y)) == "BDAB"
def test_dna() -> None:
s1 = "ACCGGTCGAGTGCGCGGAAGCCGGCCGAA"
s2 = "GTCGTTCGGAATGCCGTTGCTCTGTAAA"
assert "".join(hunt_szymanski(s1, s2)) == "GTCGTCGGAAGCCGGCCGAA"
def test_edges() -> None:
assert hunt_szymanski("", "") == []
assert hunt_szymanski("A", "") == []
assert hunt_szymanski("ABC", "DEF") == []
assert hunt_szymanski("AAAA", "AA") == ["A", "A"]
assert hunt_szymanski("ABCD", "ABCD") == list("ABCD")
assert hunt_szymanski([1, 2, 3], [3, 2, 1]) in ([1], [2], [3])
def test_random_against_table() -> None:
rng = random.Random(1977)
for _ in range(600):
x = "".join(rng.choice("ABCD"[: rng.randint(1, 4)]) for _ in range(rng.randrange(15)))
y = "".join(rng.choice("ABCD") for _ in range(rng.randrange(15)))
z = hunt_szymanski(x, y)
assert len(z) == table_lcs_length(x, y)
assert is_subsequence(z, list(x)) and is_subsequence(z, list(y))
What the next level needs
- the tests did not run in this build