Skip to content

Commit 0622991

Browse files
mitsenkovmitsenkov
committed
Handle UniProt isoforms from ALTERNATIVE PRODUCTS CC blocks
1 parent 793ddf7 commit 0622991

4 files changed

Lines changed: 361 additions & 14 deletions

File tree

‎tests/test_uniprot_extract_subset.py‎

Lines changed: 41 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
from uniprot.scripts.extract_subset import parse_de_sections
1+
from uniprot.scripts.extract_subset import build_entry_payload, parse_alt_products, parse_de_sections
22

33

44
def test_parse_de_sections_includes_subname_full() -> None:
@@ -40,3 +40,43 @@ def test_parse_de_sections_includes_alt_special_name_fields() -> None:
4040

4141
assert full_names == ["Api m 1", "Etanercept", "CD177", "Adalimumab"]
4242
assert short_names == []
43+
44+
45+
def test_parse_alt_products_splits_multiple_isoids() -> None:
46+
lines = [
47+
"CC -!- ALTERNATIVE PRODUCTS:",
48+
"CC Event=Alternative splicing; Named isoforms=3;",
49+
"CC Name=1;",
50+
"CC IsoId=P12345-1, P12345-2; Sequence=Displayed;",
51+
"CC Name=2;",
52+
"CC IsoId=P12345-3; Sequence=VSP_000001, VSP_000002;",
53+
"CC -!- SUBCELLULAR LOCATION: Cytoplasm.",
54+
]
55+
56+
isoforms = parse_alt_products(lines)
57+
58+
assert isoforms["P12345-1"] == ["Displayed"]
59+
assert isoforms["P12345-2"] == ["Displayed"]
60+
assert isoforms["P12345-3"] == ["VSP_000001", "VSP_000002"]
61+
62+
63+
def test_build_entry_payload_skips_external_isoforms() -> None:
64+
lines = [
65+
"ID TEST_HUMAN Reviewed; 5 AA.",
66+
"AC P12345;",
67+
"DE RecName: Full=Test protein;",
68+
"SQ SEQUENCE 5 AA; 500 MW; ABCDEF1234567890 CRC64;",
69+
" MTEST",
70+
"CC -!- ALTERNATIVE PRODUCTS:",
71+
"CC Event=Alternative splicing; Named isoforms=2;",
72+
"CC Name=1;",
73+
"CC IsoId=P12345-1; Sequence=Displayed;",
74+
"CC Name=2;",
75+
"CC IsoId=P12345-2; Sequence=External;",
76+
]
77+
78+
rows = build_entry_payload(lines, ["P12345"], "2025_04", reviewed=True)
79+
accessions = [row["primary_ac"] for row in rows]
80+
81+
assert "P12345-1" in accessions
82+
assert "P12345-2" not in accessions

‎uniprot/repro/isoform_audit/scripts/01_scan_sprot_isoforms_parallel.sh‎

Lines changed: 150 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,150 @@
1+
#!/usr/bin/env bash
2+
set -euo pipefail
3+
4+
SCRIPT_DIR="$(cd "$(dirname "${BASH_SOURCE[0]}")" && pwd)"
5+
REPRO_DIR="$(cd "${SCRIPT_DIR}/.." && pwd)"
6+
REPO_ROOT="$(cd "${REPRO_DIR}/../../.." && pwd)"
7+
8+
SHARDS_ROOT="/mnt/disks/toolkit-data/uniprot_shards/2025_04/sprot"
9+
JOBS=40
10+
RELEASE="2025_04"
11+
OUT_DIR="${REPRO_DIR}/outputs"
12+
LOG_DIR="${REPRO_DIR}/logs"
13+
FOCUS_ACCESSIONS="P42167 O43687 Q61029"
14+
15+
if [[ -x "${REPO_ROOT}/.venv/bin/python" ]]; then
16+
PYTHON_CMD="${REPO_ROOT}/.venv/bin/python"
17+
else
18+
PYTHON_CMD="python3"
19+
fi
20+
21+
usage() {
22+
cat <<EOF
23+
Usage: $(basename "$0") [options]
24+
25+
Parallel isoform audit for Swiss-Prot shards.
26+
27+
Options:
28+
--shards-root PATH Swiss-Prot shards dir (default: ${SHARDS_ROOT})
29+
-j, --jobs N Max parallel shard jobs (default: ${JOBS})
30+
--release TAG Release label passed to parser (default: ${RELEASE})
31+
--out-dir PATH Output dir (default: ${OUT_DIR})
32+
--log-dir PATH Log dir (default: ${LOG_DIR})
33+
--focus "A B C" Space-separated primary accessions to report in detail
34+
(default: "${FOCUS_ACCESSIONS}")
35+
--python CMD Python executable (default: ${PYTHON_CMD})
36+
-h, --help Show help
37+
EOF
38+
}
39+
40+
while [[ $# -gt 0 ]]; do
41+
case "$1" in
42+
--shards-root) SHARDS_ROOT="$2"; shift 2 ;;
43+
-j|--jobs) JOBS="$2"; shift 2 ;;
44+
--release) RELEASE="$2"; shift 2 ;;
45+
--out-dir) OUT_DIR="$2"; shift 2 ;;
46+
--log-dir) LOG_DIR="$2"; shift 2 ;;
47+
--focus) FOCUS_ACCESSIONS="$2"; shift 2 ;;
48+
--python) PYTHON_CMD="$2"; shift 2 ;;
49+
-h|--help) usage; exit 0 ;;
50+
*) echo "Unknown option: $1" >&2; usage; exit 1 ;;
51+
esac
52+
done
53+
54+
mkdir -p "${OUT_DIR}/per_shard" "${LOG_DIR}"
55+
56+
TASKS_FILE="${LOG_DIR}/sprot_isoform_tasks.txt"
57+
find "${SHARDS_ROOT}" -maxdepth 1 -type f -name 'sprot-shard-*.dat.gz' | sort > "${TASKS_FILE}"
58+
59+
if [[ ! -s "${TASKS_FILE}" ]]; then
60+
echo "[error] No shard files found in ${SHARDS_ROOT}" >&2
61+
exit 1
62+
fi
63+
64+
echo "[info] Running isoform audit across $(wc -l < "${TASKS_FILE}") Swiss-Prot shards (max parallel: ${JOBS})"
65+
cat "${TASKS_FILE}" | xargs -r -P "${JOBS}" -I {} bash -c '
66+
set -euo pipefail
67+
shard_file="{}"
68+
shard_name="$(basename "${shard_file}" .dat.gz)"
69+
out_json="'"${OUT_DIR}"'/per_shard/${shard_name}.json"
70+
log_file="'"${LOG_DIR}"'/${shard_name}.log"
71+
echo "[info] ${shard_name}" > "${log_file}"
72+
OMP_NUM_THREADS=1 OPENBLAS_NUM_THREADS=1 MKL_NUM_THREADS=1 \
73+
'"${PYTHON_CMD}"' "'"${SCRIPT_DIR}"'/scan_isoform_shard.py" \
74+
--input "${shard_file}" \
75+
--output "${out_json}" \
76+
--release "'"${RELEASE}"'" \
77+
--focus '"${FOCUS_ACCESSIONS}"' >> "${log_file}" 2>&1
78+
'
79+
80+
SUMMARY_JSON="${OUT_DIR}/summary.json"
81+
SUMMARY_TSV="${OUT_DIR}/summary.tsv"
82+
83+
"${PYTHON_CMD}" - "${OUT_DIR}/per_shard" "${SUMMARY_JSON}" "${SUMMARY_TSV}" <<'PY'
84+
import glob
85+
import json
86+
import os
87+
import sys
88+
89+
per_shard_dir, out_json, out_tsv = sys.argv[1:4]
90+
files = sorted(glob.glob(os.path.join(per_shard_dir, "*.json")))
91+
if not files:
92+
raise SystemExit("No per-shard JSON outputs found.")
93+
94+
totals = {
95+
"records_scanned": 0,
96+
"records_with_alt_products": 0,
97+
"declared_isoforms": 0,
98+
"emittable_isoforms": 0,
99+
"emitted_isoforms": 0,
100+
"missing_emittable_isoforms": 0,
101+
"non_local_isoforms": 0,
102+
"displayed_not_dash1_isoforms": 0,
103+
}
104+
focus_examples = []
105+
missing_examples = []
106+
rows = []
107+
108+
for path in files:
109+
with open(path, "r", encoding="utf-8") as handle:
110+
data = json.load(handle)
111+
shard_name = os.path.basename(path).replace(".json", "")
112+
rows.append(
113+
{
114+
"shard": shard_name,
115+
"records_with_alt_products": data["records_with_alt_products"],
116+
"declared_isoforms": data["declared_isoforms"],
117+
"emittable_isoforms": data["emittable_isoforms"],
118+
"emitted_isoforms": data["emitted_isoforms"],
119+
"missing_emittable_isoforms": data["missing_emittable_isoforms"],
120+
}
121+
)
122+
for key in totals:
123+
totals[key] += data.get(key, 0)
124+
focus_examples.extend(data.get("focus_examples", []))
125+
missing_examples.extend(data.get("missing_examples", []))
126+
127+
summary = {
128+
"totals": totals,
129+
"focus_examples": focus_examples,
130+
"missing_examples_sample": missing_examples[:100],
131+
}
132+
with open(out_json, "w", encoding="utf-8") as handle:
133+
json.dump(summary, handle, indent=2)
134+
handle.write("\n")
135+
136+
with open(out_tsv, "w", encoding="utf-8") as handle:
137+
handle.write(
138+
"shard\trecords_with_alt_products\tdeclared_isoforms\temittable_isoforms\t"
139+
"emitted_isoforms\tmissing_emittable_isoforms\n"
140+
)
141+
for row in rows:
142+
handle.write(
143+
f"{row['shard']}\t{row['records_with_alt_products']}\t{row['declared_isoforms']}\t"
144+
f"{row['emittable_isoforms']}\t{row['emitted_isoforms']}\t{row['missing_emittable_isoforms']}\n"
145+
)
146+
PY
147+
148+
echo "[info] Done."
149+
echo "[info] Summary: ${SUMMARY_JSON}"
150+
echo "[info] Per-shard table: ${SUMMARY_TSV}"

‎uniprot/repro/isoform_audit/scripts/scan_isoform_shard.py‎

Lines changed: 137 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,137 @@
1+
#!/usr/bin/env python3
2+
"""
3+
Scan one Swiss-Prot shard for ALTERNATIVE PRODUCTS isoform coverage.
4+
"""
5+
6+
from __future__ import annotations
7+
8+
import argparse
9+
import json
10+
from pathlib import Path
11+
from typing import Dict, List, Sequence
12+
13+
from uniprot.scripts.extract_subset import (
14+
build_entry_payload,
15+
extract_accessions,
16+
iter_records,
17+
parse_alt_products,
18+
parse_var_seq,
19+
)
20+
21+
22+
NON_LOCAL_SEQUENCE_SOURCES = {"external", "not described"}
23+
24+
25+
def parse_args() -> argparse.Namespace:
26+
parser = argparse.ArgumentParser(description="Audit isoform parsing for one Swiss-Prot shard.")
27+
parser.add_argument("--input", required=True, type=Path, help="Swiss-Prot shard (*.dat.gz).")
28+
parser.add_argument("--output", required=True, type=Path, help="Output JSON summary path.")
29+
parser.add_argument(
30+
"--focus",
31+
nargs="*",
32+
default=[],
33+
help="Optional list of base accessions to capture detailed examples for.",
34+
)
35+
parser.add_argument("--release", default="2025_04", help="Release label for parser call.")
36+
return parser.parse_args()
37+
38+
39+
def is_emittable_isoform(tokens: Sequence[str], varseqs: Dict[str, tuple[int, int, str]]) -> bool:
40+
if not tokens:
41+
return False
42+
lowered = [token.lower() for token in tokens]
43+
if len(tokens) == 1 and lowered[0] == "displayed":
44+
return True
45+
if any(token in NON_LOCAL_SEQUENCE_SOURCES for token in lowered):
46+
return False
47+
vsp_tokens = [token for token in tokens if token.lower() != "displayed"]
48+
return bool(vsp_tokens) and all(token in varseqs for token in vsp_tokens)
49+
50+
51+
def main() -> int:
52+
args = parse_args()
53+
focus_accessions = set(args.focus)
54+
55+
stats = {
56+
"input_file": str(args.input),
57+
"records_scanned": 0,
58+
"records_with_alt_products": 0,
59+
"declared_isoforms": 0,
60+
"emittable_isoforms": 0,
61+
"emitted_isoforms": 0,
62+
"missing_emittable_isoforms": 0,
63+
"non_local_isoforms": 0,
64+
"displayed_not_dash1_isoforms": 0,
65+
"focus_examples": [],
66+
}
67+
missing_examples: List[Dict[str, object]] = []
68+
69+
for lines in iter_records(args.input):
70+
stats["records_scanned"] += 1
71+
accessions = extract_accessions(lines)
72+
if not accessions:
73+
continue
74+
primary_ac = accessions[0]
75+
isoform_map = parse_alt_products(lines)
76+
if not isoform_map:
77+
continue
78+
79+
stats["records_with_alt_products"] += 1
80+
stats["declared_isoforms"] += len(isoform_map)
81+
varseqs = parse_var_seq(lines)
82+
payload_rows = build_entry_payload(lines, accessions, args.release, reviewed=True)
83+
emitted_isoforms = {row["primary_ac"] for row in payload_rows if row.get("is_isoform")}
84+
stats["emitted_isoforms"] += len(emitted_isoforms)
85+
86+
local_focus_details = {
87+
"primary_ac": primary_ac,
88+
"isoforms": [],
89+
}
90+
has_focus = primary_ac in focus_accessions
91+
92+
for isoform_id, tokens in isoform_map.items():
93+
lowered = [token.lower() for token in tokens]
94+
emittable = is_emittable_isoform(tokens, varseqs)
95+
emitted = isoform_id in emitted_isoforms
96+
if emittable:
97+
stats["emittable_isoforms"] += 1
98+
if not emitted:
99+
stats["missing_emittable_isoforms"] += 1
100+
if len(missing_examples) < 50:
101+
missing_examples.append(
102+
{
103+
"primary_ac": primary_ac,
104+
"isoform_id": isoform_id,
105+
"tokens": tokens,
106+
}
107+
)
108+
elif any(token in NON_LOCAL_SEQUENCE_SOURCES for token in lowered):
109+
stats["non_local_isoforms"] += 1
110+
111+
if "displayed" in lowered and not isoform_id.endswith("-1"):
112+
stats["displayed_not_dash1_isoforms"] += 1
113+
114+
if has_focus:
115+
local_focus_details["isoforms"].append(
116+
{
117+
"isoform_id": isoform_id,
118+
"tokens": tokens,
119+
"emittable": emittable,
120+
"emitted": emitted,
121+
}
122+
)
123+
124+
if has_focus:
125+
stats["focus_examples"].append(local_focus_details)
126+
127+
stats["missing_examples"] = missing_examples
128+
129+
args.output.parent.mkdir(parents=True, exist_ok=True)
130+
with args.output.open("w", encoding="utf-8") as handle:
131+
json.dump(stats, handle, indent=2)
132+
handle.write("\n")
133+
return 0
134+
135+
136+
if __name__ == "__main__":
137+
raise SystemExit(main())

0 commit comments

Comments
 (0)