Files

319 lines
8.8 KiB
Markdown
Raw Permalink Normal View History

2026-08-03 18:42:47 +08:00
# 多物种基因家族分析工作流设计
## 核心概念
### 物种角色定义
| 角色 | 说明 | 示例 |
|:---|:---|:---|
| **目标物种** (Target) | 要分析基因家族的物种 | 苦荞 (Ft) |
| **参考物种** (Reference) | 用于辅助鉴定的已知物种 | 拟南芥 (Ath)、水稻 (Os)、番茄 (Sl) |
### 文件命名规范
```
# 目标物种文件(必需)
{prefix}.pep.fasta # 目标物种蛋白序列
{prefix}.gff3 # 目标物种注释文件
{prefix}.genome.fasta # 目标物种基因组(可选)
# 参考物种文件(可选,可多个)
{ref_species}_ref.pep.fasta # 参考物种蛋白(如 Ath_ref.pep.fasta)
{ref_species}_ref.gff3 # 参考物种注释(如 Ath_ref.gff3)
# 原始数据库文件(自动识别)
Athaliana_167_TAIR10.protein_primaryTranscriptOnly.fa.gz # 拟南芥JGI格式
Osativa_323_v7.protein_primaryTranscriptOnly.fa.gz # 水稻JGI格式
```
## 工作流程
### 阶段1: 数据准备 (genome_prepare)
**输入**: 原始数据库文件(JGI/Ensembl/NCBI格式)
**处理**:
1. 识别目标物种文件(匹配prefix)
2. 识别参考物种文件(其他文件)
3. 分别处理每个物种
4. 生成标准化输出
**输出**:
```
input/{job_id}/
├── {prefix}.pep.fasta # 目标物种蛋白
├── {prefix}.gff3 # 目标物种GFF
├── {prefix}.genome.fasta # 目标物种基因组
├── Ath_ref.pep.fasta # 参考物种1蛋白
├── Ath_ref.gff3 # 参考物种1 GFF
├── Os_ref.pep.fasta # 参考物种2蛋白
└── ...
```
### 阶段2: 基因家族鉴定 (identify)
**输入**:
- 目标物种蛋白: `{prefix}.pep.fasta`
- 参考物种蛋白: `*_ref.pep.fasta`(多个)
- HMM文件或Pfam ID
**处理**:
1. HMM搜索目标物种蛋白组
2. 用参考物种蛋白构建BLAST数据库
3. 双向BLAST验证
4. 整合结果
**输出**:
```
final/{prefix}.final_members.pep.fa # 目标物种基因家族成员
final/{prefix}.identify.manifest.json # 鉴定结果清单
```
### 阶段3: 进化分析 (phylogeny)
**输入**:
- 目标物种基因家族成员
- 参考物种同源基因
**处理**:
1. 多序列比对(目标+参考物种)
2. 构建系统发育树
3. 亚家族分类
**输出**:
```
final/{prefix}.tree.nwk # 系统发育树
final/{prefix}.geneID.subfamily # 亚家族分类
```
## 关键改进点
### 1. 自动区分目标物种和参考物种
```python
class SpeciesClassifier:
"""自动分类目标物种和参考物种"""
def classify(self, files: List[Path], prefix: str) -> Tuple[SpeciesFiles, List[SpeciesFiles]]:
"""
返回: (目标物种, 参考物种列表)
"""
target = None
references = []
for file in files:
# 检查是否是目标物种
if self.is_target_species(file, prefix):
target = self.parse_target_file(file)
else:
# 参考物种
ref = self.parse_reference_file(file)
references.append(ref)
return target, references
def is_target_species(self, file: Path, prefix: str) -> bool:
"""判断文件是否属于目标物种"""
# 策略1: 文件名匹配prefix
if file.name.startswith(prefix.split('.')[0]):
return True
# 策略2: 用户明确标记
if '.target.' in file.name or '_target.' in file.name:
return True
# 策略3: 只有一个主要物种文件
return False
```
### 2. 多参考物种支持
```python
class MultiReferenceBlast:
"""多参考物种BLAST验证"""
def build_combined_db(self, reference_files: List[Path]) -> Path:
"""合并多个参考物种构建BLAST数据库"""
combined_fasta = self.work_dir / "combined_ref.pep.fa"
with open(combined_fasta, 'w') as out:
for ref_file in reference_files:
species_code = self.extract_species_code(ref_file)
with open(ref_file) as f:
for line in f:
if line.startswith('>'):
# 添加物种标识到序列ID
out.write(f">{species_code}|{line[1:]}")
else:
out.write(line)
return combined_fasta
def validate_hits(self, query_hits: List[BlastHit]) -> Dict[str, Any]:
"""验证hit来自哪个参考物种"""
results = {}
for hit in query_hits:
species = hit.target_id.split('|')[0]
if species not in results:
results[species] = []
results[species].append(hit)
return results
```
### 3. 智能物种识别
```python
class SpeciesRecognizer:
"""智能识别文件属于哪个物种"""
SPECIES_SIGNATURES = {
'Ath': {
'patterns': ['Athaliana', 'TAIR', 'AT1G', 'AT2G'],
'name': 'Arabidopsis thaliana'
},
'Os': {
'patterns': ['Osativa', 'Os01g', 'Os02g'],
'name': 'Oryza sativa'
},
'Ft': {
'patterns': ['Fagopyrum', 'FtPinG', 'tataricum'],
'name': 'Fagopyrum tataricum'
},
'Sl': {
'patterns': ['Solanum', 'Soly', 'Solyc'],
'name': 'Solanum lycopersicum'
},
'Zm': {
'patterns': ['Zea', 'Zm00001'],
'name': 'Zea mays'
}
}
def recognize_from_content(self, fasta_file: Path) -> str:
"""从FASTA文件内容识别物种"""
with open(fasta_file) as f:
# 读取前100个序列ID
ids = []
for line in f:
if line.startswith('>'):
ids.append(line[1:].split()[0])
if len(ids) >= 100:
break
# 匹配物种特征
for species_code, info in self.SPECIES_SIGNATURES.items():
match_count = sum(1 for id in ids if any(p in id for p in info['patterns']))
if match_count > len(ids) * 0.5: # 超过50%匹配
return species_code
return "UNKNOWN"
```
## 实现建议
### 方案A: 文件名约定(简单)
用户按约定命名文件:
```
# 目标物种
Ft.pep.fasta # 必需
Ft.gff3 # 可选
# 参考物种
Ath_ref.pep.fasta # 拟南芥参考
Os_ref.pep.fasta # 水稻参考
Sl_ref.pep.fasta # 番茄参考
```
### 方案B: 配置文件(灵活)
```yaml
# species_config.yaml
target:
code: "Ft"
name: "Fagopyrum tataricum"
files:
pep: "Ft.pep.fasta"
gff: "Ft.gff3"
references:
- code: "Ath"
name: "Arabidopsis thaliana"
files:
pep: "Athaliana_167_TAIR10.protein_primaryTranscriptOnly.fa.gz"
gff: "Athaliana_167_TAIR10.gene_exons.gff3.gz"
- code: "Os"
name: "Oryza sativa"
files:
pep: "Osativa_323_v7.protein_primaryTranscriptOnly.fa.gz"
- code: "Sl"
name: "Solanum lycopersicum"
files:
pep: "Slycopersicum_...fa.gz"
```
### 方案C: 自动识别(智能)
系统自动识别:
1. 扫描所有文件
2. 从内容识别物种
3. 用户确认或修正
## 推荐实现
**短期**(本周):方案A - 文件名约定
- 简单,用户容易理解
- 改动最小
**中期**(本月):方案B - 配置文件
- 更灵活,支持复杂场景
- 可扩展性好
**长期**(未来):方案C - 自动识别+LLM
- 完全自动化
- 智能推荐参考物种
## 立即实施的修改
### 1. 修改 `genome_prepare`
```python
def main():
parser.add_argument("--target-prefix", required=True, help="目标物种前缀")
parser.add_argument("--reference-species", nargs='+', help="参考物种代码列表")
# 处理逻辑
target_files = find_target_files(input_dir, args.target_prefix)
reference_files = find_reference_files(input_dir, exclude=target_files)
# 分别处理
process_target(target_files)
for ref in reference_files:
process_reference(ref)
```
### 2. 修改 `identify`
```python
def run_identify():
# 加载目标物种
target_pep = load_target_protein(prefix)
# 加载参考物种
reference_peps = load_reference_proteins(workspace)
# HMM搜索目标物种
hmm_results = run_hmmsearch(target_pep)
# BLAST验证(使用所有参考物种)
combined_db = build_combined_blast_db(reference_peps)
blast_results = run_blastp(target_pep, combined_db)
# 整合结果
final_members = integrate_results(hmm_results, blast_results)
```
这样设计是否满足你的需求?需要我立即实施哪个方案?