Skip to content

Commit 0c262a6

Browse files
committed
UPD
1 parent 5390d6b commit 0c262a6

1 file changed

Lines changed: 139 additions & 0 deletions

File tree

scripts/filter_gbk_by_product.py

Lines changed: 139 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,139 @@
1+
"""Filter GenBank files by antiSMASH cand_cluster/protocluster product and copy matches."""
2+
3+
from __future__ import annotations
4+
5+
import argparse
6+
import os
7+
import re
8+
import shutil
9+
import sys
10+
from concurrent.futures import ProcessPoolExecutor
11+
from pathlib import Path
12+
from typing import Iterable, Sequence
13+
14+
from tqdm import tqdm
15+
16+
17+
def _iter_gbk_files(root: Path, exts: Sequence[str]) -> Iterable[Path]:
18+
for ext in exts:
19+
yield from root.rglob(f"*{ext}")
20+
21+
22+
def _file_matches(
23+
gbk_path: str,
24+
products: Sequence[str],
25+
) -> dict | None:
26+
"""
27+
Return match info if the GBK contains cand_cluster/protocluster with desired product.
28+
Uses lightweight text search to avoid full parsing overhead.
29+
"""
30+
path = Path(gbk_path)
31+
try:
32+
text = path.read_text(errors="ignore")
33+
except Exception as exc: # noqa: BLE001
34+
return {"error": f"read failed: {exc}", "file": gbk_path}
35+
36+
lower = text.lower()
37+
if "cand_cluster" not in lower and "protocluster" not in lower:
38+
return None
39+
40+
matched: set[str] = set()
41+
for prod in products:
42+
pat = re.compile(r'/product="[^"]*' + re.escape(prod.lower()) + r'[^"]*"', re.IGNORECASE)
43+
if pat.search(text):
44+
matched.add(prod)
45+
46+
if not matched:
47+
return None
48+
49+
return {
50+
"file": path.name,
51+
"path": str(path),
52+
"matched_products": sorted(matched),
53+
}
54+
55+
56+
def cli() -> argparse.Namespace:
57+
parser = argparse.ArgumentParser(description=__doc__)
58+
parser.add_argument("input_dir", type=Path, help="Directory containing GenBank files")
59+
parser.add_argument("output_dir", type=Path, help="Directory to copy matching GBKs into")
60+
parser.add_argument(
61+
"-p",
62+
"--products",
63+
nargs="+",
64+
required=True,
65+
help='Product substrings to match (e.g., "NRPS", "T1PKS")',
66+
)
67+
parser.add_argument(
68+
"-w",
69+
"--workers",
70+
type=int,
71+
default=os.cpu_count() or 1,
72+
help="Number of worker processes",
73+
)
74+
parser.add_argument(
75+
"--exts",
76+
nargs="+",
77+
default=[".gbk", ".gbff", ".gb"],
78+
help="File extensions to include",
79+
)
80+
parser.add_argument(
81+
"--chunksize",
82+
type=int,
83+
default=16,
84+
help="Chunksize for process pool mapping",
85+
)
86+
return parser.parse_args()
87+
88+
89+
def main() -> None:
90+
args = cli()
91+
input_dir: Path = args.input_dir
92+
output_dir: Path = args.output_dir
93+
products = [p.lower() for p in args.products]
94+
95+
if not input_dir.exists():
96+
raise SystemExit(f"Input directory does not exist: {input_dir}")
97+
98+
files = list(_iter_gbk_files(input_dir, args.exts))
99+
if not files:
100+
raise SystemExit("No GenBank files found.")
101+
102+
output_dir.mkdir(parents=True, exist_ok=True)
103+
104+
total = len(files)
105+
copied = 0
106+
errors = 0
107+
108+
with ProcessPoolExecutor(max_workers=args.workers) as executor:
109+
with tqdm(total=total, desc="Scanning GBKs", unit="file") as pbar:
110+
for result in executor.map(
111+
_file_matches,
112+
(str(p) for p in files),
113+
[products] * len(files),
114+
chunksize=args.chunksize,
115+
):
116+
pbar.update(1)
117+
118+
if not result:
119+
continue
120+
121+
if isinstance(result, dict) and "error" in result:
122+
errors += 1
123+
pbar.set_postfix(copied=copied, errors=errors)
124+
continue
125+
126+
src_path = Path(result["path"])
127+
dst_path = output_dir / src_path.name
128+
try:
129+
shutil.copy2(src_path, dst_path)
130+
copied += 1
131+
except Exception as exc: # noqa: BLE001
132+
errors += 1
133+
result = {"error": f"copy failed: {exc}", "file": str(src_path)}
134+
135+
pbar.set_postfix(copied=copied, errors=errors)
136+
137+
138+
if __name__ == "__main__":
139+
main()

0 commit comments

Comments
 (0)