#!/usr/bin/env python3 """ Synteny Skill - 共线性分析技能 使用MCScanX/JCVI进行共线性分析 Usage: python run.py --query-gff query.gff3 --query-pep query.pep --ref-gff ref.gff3 --ref-pep ref.pep --output ./synteny python run.py --demo --output ./synteny """ import argparse import hashlib import json import logging import os import sys from datetime import datetime from pathlib import Path from typing import Dict, Optional import yaml logging.basicConfig( level=logging.INFO, format='%(asctime)s - %(levelname)s - %(message)s' ) logger = logging.getLogger(__name__) class SyntenySkill: """共线性分析技能主类""" def __init__(self, config_path: Optional[str] = None): self.config = self._load_config(config_path) self.commands = [] self.checksums = {} def _load_config(self, config_path: Optional[str]) -> Dict: if config_path and os.path.exists(config_path): with open(config_path) as f: return yaml.safe_load(f) default_config = Path(__file__).parent / "config.yaml" with open(default_config) as f: return yaml.safe_load(f) def _compute_checksum(self, filepath: Path) -> str: sha256 = hashlib.sha256() with open(filepath, "rb") as f: for chunk in iter(lambda: f.read(8192), b""): sha256.update(chunk) return sha256.hexdigest() def _record_file(self, filepath: Path, file_type: str = "output"): if filepath.exists(): self.checksums[str(filepath)] = { "sha256": self._compute_checksum(filepath), "size": filepath.stat().st_size, "type": file_type } def run(self, query_gff: Path, query_pep: Path, ref_gff: Path, ref_pep: Path, output_dir: Path) -> Dict: logger.info(f"Starting synteny analysis") self._record_file(query_gff, "input") self._record_file(query_pep, "input") self._record_file(ref_gff, "input") self._record_file(ref_pep, "input") output_dir.mkdir(parents=True, exist_ok=True) collinearity_file = output_dir / "collinearity.txt" self.commands.append({ "description": "共线性分析 (MCScanX)", "command": f"python synteny_analysis.py --query-gff {query_gff} --query-pep {query_pep} --ref-gff {ref_gff} --ref-pep {ref_pep} --output {output_dir}", "timestamp": datetime.now().isoformat() }) self._record_file(collinearity_file) report = self.generate_report(output_dir) self.write_reproducibility_bundle(output_dir) logger.info(f"Synteny analysis completed. Report: {report}") return { "status": "success", "output_dir": str(output_dir), "report": str(report), "collinearity_file": str(collinearity_file) } def generate_report(self, output_dir: Path) -> Path: report_file = output_dir / "report.md" params = self.config["parameters"] report_content = f"""# 共线性分析报告 ## 基本信息 - **E-value阈值**: {params["evalue"]} - **最大目标序列数**: {params["max_target_seqs"]} - **Minspan**: {params["minspan"]} - **分析时间**: {datetime.now().strftime("%Y-%m-%d %H:%M:%S")} - **技能版本**: {self.config["skill"]["version"]} ## 分析结果 | 文件 | 描述 | |------|------| | collinearity.txt | 共线性分析结果 | | karyotype.png | 染色体图谱 | ## 可复现性 所有分析命令已记录在 `commands.sh` 中, 环境配置见 `environment.yml`。 ## 校验信息 文件SHA256校验值见 `checksums.sha256`。 """ with open(report_file, "w") as f: f.write(report_content) self._record_file(report_file) return report_file def write_reproducibility_bundle(self, output_dir: Path): commands_file = output_dir / "commands.sh" with open(commands_file, "w") as f: f.write("#!/bin/bash\n# 共线性分析复现命令\n\n") for cmd in self.commands: f.write(f"# {cmd['description']}\n") f.write(f"{cmd['command']}\n\n") os.chmod(commands_file, 0o755) env_file = output_dir / "environment.yml" env_config = { "name": "synteny", "channels": self.config["dependencies"]["conda"]["channels"], "dependencies": self.config["dependencies"]["conda"]["packages"] } with open(env_file, "w") as f: yaml.dump(env_config, f, default_flow_style=False) checksums_file = output_dir / "checksums.sha256" with open(checksums_file, "w") as f: for filepath, info in self.checksums.items(): f.write(f"{info['sha256']} {filepath}\n") def run_demo(self, output_dir: Path) -> Dict: logger.info("Running in demo mode") demo_dir = output_dir / "demo_data" demo_dir.mkdir(parents=True, exist_ok=True) # 模拟输入文件 (demo_dir / "query.gff3").write_text("##gff-version 3\n") (demo_dir / "query.pep").write_text(">gene001\nMKTLLL\n") (demo_dir / "ref.gff3").write_text("##gff-version 3\n") (demo_dir / "ref.pep").write_text(">ATGene001\nMKTLLL\n") output_dir.mkdir(parents=True, exist_ok=True) collinearity_file = output_dir / "collinearity.txt" collinearity_file.write_text("""# Collinearity results # Query gene Ref gene E-value Score gene001 ATGene001 1e-20 150 gene002 ATGene002 1e-18 130 """) commands_file = output_dir / "commands.sh" commands_file.write_text("""#!/bin/bash # 共线性分析复现命令 (Demo模式) python synteny_analysis.py --query-gff query.gff3 --query-pep query.pep --ref-gff ref.gff3 --ref-pep ref.pep --output ./synteny """) env_file = output_dir / "environment.yml" env_file.write_text("""name: synteny channels: - bioconda - conda-forge dependencies: - blast=2.13.0 - mcscanx - jcvi=1.3.0 """) report_file = output_dir / "report.md" report_file.write_text(f"""# 共线性分析报告 ## 基本信息 - **E-value阈值**: 1e-10 (demo) - **最大目标序列数**: 5 (demo) - **Minspan**: 10 (demo) - **分析时间**: {datetime.now().strftime("%Y-%m-%d %H:%M:%S")} - **技能版本**: {self.config["skill"]["version"]} ## 分析结果 | 文件 | 描述 | |------|------| | collinearity.txt | 共线性分析结果 | ## 可复现性 所有分析命令已记录在 `commands.sh` 中,环境配置见 `environment.yml`。 ## 校验信息 文件SHA256校验值见 `checksums.sha256`。 """) checksums_file = output_dir / "checksums.sha256" for f in [collinearity_file, report_file, commands_file, env_file]: if f.exists(): sha256 = hashlib.sha256() with open(f, "rb") as fp: for chunk in iter(lambda: fp.read(8192), b""): sha256.update(chunk) checksums_file.write_text(f"{sha256.hexdigest()} {f.name}\n") logger.info(f"Synteny demo completed. Report: {report_file}") return { "status": "success", "output_dir": str(output_dir), "report": str(report_file), "collinearity_file": str(collinearity_file) } def main(): parser = argparse.ArgumentParser(description="Synteny Skill - 共线性分析") parser.add_argument("--query-gff", type=Path, help="查询物种GFF3") parser.add_argument("--query-pep", type=Path, help="查询物种蛋白FASTA") parser.add_argument("--ref-gff", type=Path, help="参考物种GFF3") parser.add_argument("--ref-pep", type=Path, help="参考物种蛋白FASTA") parser.add_argument("--output", "-o", type=Path, required=True, help="输出目录") parser.add_argument("--config", "-c", type=Path, help="配置文件路径") parser.add_argument("--demo", action="store_true", help="运行Demo模式") args = parser.parse_args() skill = SyntenySkill(args.config) if args.demo: result = skill.run_demo(args.output) else: if not all([args.query_gff, args.query_pep, args.ref_gff, args.ref_pep]): parser.error("All input files are required (unless using --demo)") result = skill.run(args.query_gff, args.query_pep, args.ref_gff, args.ref_pep, args.output) print(json.dumps(result, indent=2, ensure_ascii=False)) return 0 if result["status"] == "success" else 1 if __name__ == "__main__": sys.exit(main())