rasusa 5.1.0

Randomly subsample reads or alignments
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105

#!/usr/bin/env python3
"""Helpers for benches/update_readme.sh: render benchmark blocks and splice them into
README.md between the `<!-- BENCH:<name>:START/END -->` markers.
"""
import json
import re
import sys


def _relative_str(mean, stddev, fastest_mean, fastest_stddev, with_stddev=True):
    ratio = mean / fastest_mean
    if with_stddev:
        rel_err = ratio * ((stddev / mean) ** 2 + (fastest_stddev / fastest_mean) ** 2) ** 0.5
        return f"{ratio:.2f} ± {rel_err:.2f}"
    return f"{ratio:.2f}"


def _load_results(json_path):
    with open(json_path) as f:
        return json.load(f)["results"]


# Matches the rasusa invocation specifically (executable named "rasusa" followed by a
# subcommand), not just any command whose *arguments* happen to contain the substring
# "rasusa" - which they always do in CI, since the repo (and so the checkout path,
# e.g. /home/runner/work/rasusa/rasusa/...) is named "rasusa" too.
_RASUSA_CMD_RE = re.compile(r"(?:^|/)rasusa (reads|aln)\b")


def _is_rasusa_cmd(command):
    return _RASUSA_CMD_RE.search(command) is not None


def _read_table(md_path):
    with open(md_path) as f:
        return f.read().rstrip("\n")


def cmd_single_block(args):
    md_path, json_path = args
    results = _load_results(json_path)
    rasusa = next(r for r in results if _is_rasusa_cmd(r["command"]))
    other = next(r for r in results if r is not rasusa)

    table = _read_table(md_path)
    relative = _relative_str(other["mean"], other["stddev"], rasusa["mean"], rasusa["stddev"])
    print(table)
    print()
    print(f"**Summary**: `rasusa` ran {relative} times faster than `filtlong`.")


def cmd_paired_block(args):
    md_path, json_path = args
    results = _load_results(json_path)
    rasusa = next(r for r in results if _is_rasusa_cmd(r["command"]))
    others = [r for r in results if r is not rasusa]
    # Order: 1-pass seqtk, 2-pass seqtk, matching the hyperfine invocation order.
    ratio1 = _relative_str(others[0]["mean"], others[0]["stddev"], rasusa["mean"], rasusa["stddev"], with_stddev=False)
    ratio2 = _relative_str(others[1]["mean"], others[1]["stddev"], rasusa["mean"], rasusa["stddev"], with_stddev=False)

    table = _read_table(md_path)
    print(table)
    print()
    print(
        f"**Summary**: `rasusa reads` ran {ratio1} times faster than `seqtk` (1-pass) "
        f"and {ratio2} times faster than `seqtk` (2-pass)"
    )


def cmd_splice(args):
    readme_path, marker, content_path = args
    with open(readme_path) as f:
        readme = f.read()
    with open(content_path) as f:
        content = f.read().rstrip("\n")

    start = f"<!-- BENCH:{marker}:START -->"
    end = f"<!-- BENCH:{marker}:END -->"
    pattern = re.compile(re.escape(start) + r".*?" + re.escape(end), re.DOTALL)
    if not pattern.search(readme):
        print(f"error: markers {start!r}/{end!r} not found in {readme_path}", file=sys.stderr)
        sys.exit(1)

    replacement = f"{start}\n{content}\n{end}"
    readme = pattern.sub(lambda _: replacement, readme, count=1)
    with open(readme_path, "w") as f:
        f.write(readme)


COMMANDS = {
    "single-block": cmd_single_block,
    "paired-block": cmd_paired_block,
    "splice": cmd_splice,
}


def main():
    if len(sys.argv) < 2 or sys.argv[1] not in COMMANDS:
        print(f"usage: {sys.argv[0]} <{'|'.join(COMMANDS)}> [args...]", file=sys.stderr)
        sys.exit(2)
    COMMANDS[sys.argv[1]](sys.argv[2:])


if __name__ == "__main__":
    main()