基因组注释计划蓝图

gene-anno 工程蓝图

适用场景:任意物种的基因组结构 + 功能注释工程,用于从基因组序列到完整发布集的系统化加工。


目录

  1. 工程定位与总体流程
  2. 顶层目录结构
  3. 配置文件结构说明(config.yaml 占位)
  4. 必备软件与数据库
  5. 每步脚本说明(01–12 + 可选模块)
    • 5.1 公共工具模块
    • 5.2–5.13:01–12 主线
    • 5.14:03b 长读转录本证据模块(可选)
    • 5.15:10b ncRNA 注释模块(可选)
    • 5.16:12b 交叉流水线与近缘物种对照 QC(高级模块)
  6. 功能注释整合机制(本地 vs 线上)
  7. 日志规范(文件 + 终端流式输出)
  8. QC 与交付物
  9. 方法学模板(占位说明)
  10. 施工注意事项(硬件与数据组织)

1. 工程定位与总体流程

1.1 项目定位

  • 项目名称:gene-anno
  • 核心目标:为任意一个物种,提供一套可复用、工程化的基因组注释流程,包括:
    • 重复序列注释
    • RNA-seq 支持的基因结构预测
    • 可选:长读(Iso-Seq / Nanopore)转录本证据整合
    • 基因集质量控制与过滤(含 TE 相关过滤、证据打分、AED 分布、白名单回捞与假基因分流)
    • 同源、结构域、KEGG、eggNOG 等多源功能注释
    • 可选:ncRNA 注释(tRNA、rRNA、snRNA、miRNA 等)
    • 标准发布集构建与 QC 汇总
    • 可选:与其它注释流水线结果、近缘物种注释结果进行对照 QC

这一套流程强调配置驱动、目录清晰、表头契约,在不同物种、不同项目中可以重复使用。

1.2 在整体科研工程中的角色

对同一物种:

  • gene-anno

    负责产生统一的注释发布集:

    • genome.fa
    • genes.gff3
    • cds.fa
    • pep.fa
    • tx2gene.tsv
    • functional.tsv
    • TS.ncrna.gff3(如启用 ncRNA 模块)
    • 各类 QC 摘要与比较 QC 报告
  • phylo

    • 从 gene-anno 的 genome.fa + genes.gff3 中,提取规范化 CDS/pep,用于跨物种正交基因构建与系统发育分析。
  • aphylo

    • 利用 gene-anno 产生的蛋白和功能注释,进行 PSG、CAFE、MCMCTree 等进化分析。
  • rna-seq

    • 使用 gene-anno 的
      genome.fa + genes.gff3 + tx2gene.tsv + functional.tsv
      进行表达定量与富集分析。

核心思想
对每个物种,只有一套权威注释来自 gene-anno,所有下游工程都以此为统一入口,避免多版本混乱。

1.3 设计原则

  1. 可复用
    • 换物种时,只需替换 config/config.yaml 中的物种信息与 data/ 下的原始数据。
    • 脚本逻辑保持稳定,由配置和表头驱动行为。
  2. 契约清晰
    • 每一步的输入输出路径、表头格式,在蓝图中写明。
    • 脚本之间只通过这些“契约文件”通信。
  3. 可追溯
    • 每一步都有日志文件,并在终端输出关键信息。
    • 关键统计统一写入 TSV,方便查阅与复现。
  4. 可扩展
    • 在不破坏 01–12 主线的前提下,引入 03b、10b、12b 等模块;
    • 模块通过配置开关控制,按需启用。
  5. 多层 QC
    • 结构层:基因长度、外显子数、TE 重叠等;
    • 证据层:BUSCO、同源、InterPro、RNA-seq 与长读支持等,形成证据打分(E-score);
    • AED(Annotation Edit Distance)分布:量化基因模型与转录证据吻合度;
    • 高级策略:白名单回捞(SwissProt/InterPro 强支持)、假基因分流、比较 QC。

1.4 12 步流程总览 + 可选模块

主线步骤顺序:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
01_genome_qc
02_repeat_annot
03_rnaseq_qc_and_clean
03b_longread_support # 可选:长读 / Iso-Seq 证据整合
04_rnaseq_mapping
05_braker3_run
06_gene_qc_and_filter
07_functional_homology
08_functional_interpro
09_functional_kegg # 支持本地 / 线上
10_functional_eggnog # 支持本地 / 线上
10b_ncrna_annot # 可选:ncRNA 注释
11_build_publish_release
12_novel_gene # 可选:新基因候选
12b_comparative_qc # 可选:交叉流水线 + 近缘物种对照 QC

可以理解为四层结构:

  • 基础层(01–04 + 03b)
    • 基因组 QC 与重复注释;
    • RNA-seq 质控与比对;
    • 如有长读,整合转录本证据,生成 hints。
  • 基因集层(05–06)
    • 利用 masked genome + RNA-seq + 长读 + 蛋白证据进行结构预测;
    • 对基因集进行统计、TE 过滤、证据打分、AED 分析;
    • 引入白名单回捞机制与假基因分流,生成高可信基因列表。
  • 功能与发布层(07–12 + 10b)
    • 多源功能注释(同源、InterPro、KEGG、eggNOG);
    • 可选 ncRNA 注释;
    • 整合功能表,构建发布集;
    • 可选新基因候选筛选。
  • 高级 QC 层(12b)
    • 与其它流水线结果、近缘物种注释对照;
    • 汇总比较 QC 指标,辅助判断注释整体合理性。

1.5 标准发布集内容

以物种代码 TS 为例,完成 gene-anno 后至少产出:

  • 结构与序列:
    • TS.genome.fa
    • TS.genes.gff3
    • TS.cds.fa
    • TS.pep.fa
    • TS.tx2gene.tsv
  • 功能与 QC:
    • TS.functional.tsv
    • TS.busco_proteome_summary.txt
    • TS.qc_summary.tsv
    • README.annotation.txt
  • 可选扩展:
    • TS.ncrna.gff3:ncRNA 注释(10b)
    • results/12b_comparative_qc/:交叉流水线与近缘物种对照 QC 报告(12b)
    • results/06_gene_qc/aed_distribution.tsv:AED 分布表(建议)
    • results/06_gene_qc/gene.confidence.tsv:证据打分与标签表

这些文件作为 phylo / aphylo / rna-seq / NCBI 提交 / 论文附录的统一基础。


2. 顶层目录结构

2.1 目录树与各部分职责

推荐顶层结构:

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
gene-anno/
├── config/
│ └── config.yaml
├── scripts/
│ ├── common_utils.py
│ ├── 01_genome_qc.py
│ ├── 02_repeat_annot.py
│ ├── 03_rnaseq_qc_and_clean.py
│ ├── 03b_longread_support.py # 可选长读转录本模块
│ ├── 04_rnaseq_mapping.py
│ ├── 05_braker3_run.py
│ ├── 06_gene_qc_and_filter.py
│ ├── 07_functional_homology.py
│ ├── 08_functional_interpro.py
│ ├── 09_functional_kegg.py
│ ├── 10_functional_eggnog.py
│ ├── 10b_ncrna_annot.py # 可选 ncRNA 注释模块
│ ├── 11_build_publish_release.py
│ ├── 12_novel_gene.py
│ └── 12b_comparative_qc.py # 可选高级 QC 模块
├── data/
│ ├── genome/
│ ├── rnaseq/
│ ├── longread/ # 可选:Iso-Seq / Nanopore 转录本
│ ├── proteins/
│ ├── annotation/ # TE / 外部 ncRNA / 外部 GFF 注释
│ └── functional_external/ # 线上 KEGG/eggNOG 结果(KAAS/eggNOG-mapper)
├── db/
│ ├── busco/
│ ├── dfam/
│ ├── rfam/
│ ├── interpro/
│ ├── nr/
│ ├── swissprot/
│ ├── kegg/
│ └── eggnog/
├── results/
│ ├── 01_genome_qc/
│ ├── 02_repeat/
│ ├── 03_rnaseq_qc/
│ ├── 03b_longread_support/
│ ├── 04_mapping/
│ ├── 05_braker/
│ ├── 06_gene_qc/
│ ├── 07_homology/
│ ├── 08_interpro/
│ ├── 09_kegg/
│ ├── 10_eggnog/
│ ├── 10b_ncrna/ # 可选:ncRNA 结果
│ ├── 11_publish/
│ ├── 12_novel_gene/
│ └── 12b_comparative_qc/ # 可选:对照 QC 结果
└── logs/
  • config/:全局配置。
  • scripts/:步骤脚本与公共模块。
  • data/:基因组、RNA-seq、长读、外部注释与线上功能注释结果。
  • db/:共享数据库。
  • results/:各步骤产物。
  • logs/:所有脚本运行日志。

2.2 data/ 子目录使用原则

  • data/genome/:原始与 clean 基因组。
  • data/rnaseq/:原始 / clean RNA-seq 及 meta 表。
  • data/longread/:长读原始 reads、BAM、GTF。
  • data/proteins/:近缘物种蛋白。
  • data/annotation/:外部 TE/ncRNA/GFF 注释。
  • data/functional_external/:KAAS、eggNOG-mapper 等线上结果表。

3. 配置文件结构说明(config.yaml 占位)

此处仅占位,详情见附录中完整 的config.yaml。


4. 必备软件与数据库

4.1 软件组件说明

与原蓝图一致,只补充一点约定:

  1. 质控与比对
    • fastp / FastQC
    • HISAT2 / STAR
    • samtools
    • minimap2
  2. 结构注释
    • BRAKER3(集成 GeneMark、AUGUSTUS 等)
  3. 重复与 ncRNA
    • RepeatModeler / RepeatMasker / EDTA
    • Infernal(cmscan)+ Rfam
    • tRNAscan-SE
  4. 功能注释
    • DIAMOND
    • InterProScan
    • KEGG 注释工具或对 KAAS 结果解析脚本
    • eggNOG-mapper 或对其 web 结果解析脚本
  5. 质量评估与辅助
    • BUSCO
    • bedtools / gffread / gffutils

4.2 数据库部署

类型 内容示例 路径示例
BUSCO lineage mollusca_odb10 db/busco/
TE 库 Dfam / EDTA 库 db/dfam/
ncRNA Rfam CM db/rfam/
Domain / GO InterProScan DB db/interpro/
同源注释(精) SwissProt db/swissprot/
同源注释(全) NCBI NR db/nr/
KEGG(可选) KEGG KO DB db/kegg/
eggNOG(可选) eggNOG DB db/eggnog/

5. 每步脚本说明(01–12 + 可选模块)

5.1 common_utils.py —— 公共工具模块

目标:统一基础能力:

  • 读取 config.yaml
  • 构造 logger(写日志 + 终端输出)
  • 封装外部程序调用(失败时记录命令、退出码和 stderr 摘要)
  • 统一 TSV 读写,检查必备表头
  • 建立目录、检查输入存在、写简单统计信息

日志规范约定:

  • 每脚本对应 logs/XX_stepname.log
  • 日志中至少包含:
    • config 版本摘要(如 git tag / 手动版本号)
    • 本次运行使用的关键参数(数据库、线程数等)
    • 输入文件路径与记录数
    • 核心统计(例如 BUSCO C/F/M、基因数等)
    • 错误与警告(ERROR / WARNING)

5.2 01_genome_qc.py —— 基因组质量评估

目的
评估基因组组装质量,为后续注释提供背景。

输入

  • paths.genome.clean
  • db.busco_lineage

输出

  • results/01_genome_qc/genome.busco/
  • results/01_genome_qc/genome.stats.tsv

表头约定:

genome.stats.tsv

1
2
3
metric      # 指标名称,如 total_length, n_contigs, n50, gc_content, gap_count
value # 数值
note # 说明或单位,如 "bp", "fraction"

BUSCO 输出使用官方格式,不额外规定。


5.3 02_repeat_annot.py —— 重复序列注释

目的
识别并标记基因组中的重复序列,生成 TE GFF 与 masked genome。

输入

  • paths.genome.clean
  • db.dfam 或 EDTA 内部库

输出

  • results/02_repeat/TS.repeat.lib.fa
  • results/02_repeat/TS.repeats.gff3
  • results/02_repeat/TS.genome.masked.fa
  • results/02_repeat/repeat.stats.tsv

表头约定:

repeat.stats.tsv

1
2
3
4
class        # 重复类型,如 LINE, SINE, LTR, DNA, Simple_repeat 等
n_elements # 元件数量
total_bp # 覆盖的碱基总长
genome_frac # 占基因组比例(0-1)

TS.repeats.gff3 遵循 GFF3 标准,attributes 中推荐包含:

  • ID
  • Name
  • Class or Family

5.4 03_rnaseq_qc_and_clean.py —— RNA-seq 质控与清洗

目的
对 RNA-seq 样本进行统一 QC,生成 clean reads。

输入

  • data/rnaseq/raw/*.fq.gz
  • data/rnaseq/samples.meta.tsv

样本 meta 表必需列:

1
2
3
4
5
6
7
sample_id
run_accession
library_layout # PE/SE
read_length
tissue
condition
replicate

输出

  • data/rnaseq/clean/:clean fq
  • results/03_rnaseq_qc/fastp_summary.tsv

表头约定:

fastp_summary.tsv

1
2
3
4
5
6
7
8
9
sample_id
raw_reads
clean_reads
raw_bases
clean_bases
q30_rate
gc_content
adapter_trimmed_reads
note

5.5 03b_longread_support.py —— 长读转录本证据整合(可选)

目的
将 Iso-Seq / Nanopore 转录本证据转换为 BRAKER3 可用的 hints。

输入

  • paths.longread.bampaths.longread.gtf
  • 可选:longread_transcripts.fa

输出

  • results/03b_longread_support/longread_hints.gff
  • results/03b_longread_support/longread_transcripts.stats.tsv

表头约定:

longread_transcripts.stats.tsv

1
2
3
4
5
6
7
transcript_id
source # isoseq3 / tama / stringtie2 / minimap2
n_exons
length_nt
n_supporting_reads
has_polya
support_level # high / medium / low

5.6 04_rnaseq_mapping.py —— RNA-seq 比对

目的
将 clean reads 比对到基因组,生成 BAM 与 splice hints。

输入

  • data/rnaseq/clean/
  • paths.genome.clean 或 masked genome

输出

  • results/04_mapping/rnaseq.sorted.bam(或多个 BAM)
  • results/04_mapping/hints.gff(可选)
  • results/04_mapping/mapping.stats.tsv

表头约定:

mapping.stats.tsv

1
2
3
4
5
6
7
8
sample_id
total_reads
mapped_reads
uniquely_mapped_reads
mapping_rate
unique_mapping_rate
mean_coverage
note

5.7 05_braker3_run.py —— 基因结构注释

目的
在重复屏蔽基因组上,整合 RNA-seq、长读和蛋白证据进行基因预测。

输入

  • results/02_repeat/TS.genome.masked.fa
  • results/04_mapping/rnaseq.sorted.bam + hints.gff
  • results/03b_longread_support/longread_hints.gff(如启用)
  • data/proteins/homologs.merged.pep.fa

输出

  • results/05_braker/braker.gff3
  • results/05_braker/braker.cds.fa
  • results/05_braker/braker.proteins.fa
  • results/05_braker/braker.log.summary.tsv

表头约定(log.summary)

1
2
3
metric
value
note

例如:

  • metric = n_genes, value = 32000
  • metric = n_transcripts, value = xxx 等。

5.8 06_gene_qc_and_filter.py —— 基因集 QC 与过滤

目的
对 BRaker 基因集进行结构统计、TE 过滤、证据评分和 AED 分析,形成高可信基因列表、假基因列表与低可信基因列表。

输入

  • braker.gff3, braker.cds.fa, braker.proteins.fa
  • results/02_repeat/TS.repeats.gff3
  • results/04_mapping/rnaseq.sorted.bam(用于 AED)
  • 可选:长读 BAM/转录本 GTF(用于 AED)
  • results/07_homology / results/08_interpro / BUSCO 结果(供白名单回捞与证据评分使用)

注:证据打分 / AED 统计可分多次运行,蓝图允许脚本内部串联。

输出(核心):

  • results/06_gene_qc/gene.stats.tsv
  • results/06_gene_qc/gene.confidence.tsv
  • results/06_gene_qc/gene.filtered_gene_ids.txt(高可信基因 ID)
  • results/06_gene_qc/pseudogene.gff3(假基因)
  • results/06_gene_qc/aed_distribution.tsv
  • results/06_gene_qc/proteome.busco/

表头约定:

gene.stats.tsv

1
2
3
4
5
6
7
8
9
gene_id
transcript_id
cds_len_nt
pep_len_aa
n_exon
gc_cds
te_overlap_fraction
is_single_exon
busco_status # C/F/M/NA 或具体标签

gene.confidence.tsv(证据层 + 策略层):

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
gene_id
cds_len_nt
pep_len_aa
te_overlap_fraction
has_rnaseq_support # true/false
has_longread_support # true/false
has_swissprot_hit # true/false
has_interpro_domain # true/false
has_eggnog_og # true/false
busco_status # C/F/M/NA
evidence_score # E-score(整数或浮点)
aed # 0-1
is_te_like # true/false
is_short_cds # true/false
has_frame_shift_or_stop # true/false
is_pseudogene # true/false(结构损坏 + 同源支持)
is_rescued_by_homology # true/false(白名单回捞)
keep_flag # keep / drop
confidence_level # high / medium / low

aed_distribution.tsv

1
2
3
aed_bin       # 如 [0.0-0.1), [0.1-0.2) 等
n_genes
fraction

过滤与策略要点(概念):

  • 结构过滤:
    • 短 CDS(min_cds_len_ntmin_protein_len_aa
    • 高 TE 重叠(te_overlap_fraction
    • 过多 N 或 frameshift
  • 证据评分(E-score):
    • 来自 RNA-seq、长读、SwissProt、InterPro、eggNOG、BUSCO 等;
    • 脚本中定义权重,不在蓝图展开;
  • AED 分布:
    • 输出 aed_distribution.tsv 用于整体查看;
  • 白名单回捞:
    • 即使基因短或高 TE 重叠,只要 has_swissprot_hithas_interpro_domain 等强证据,可 is_rescued_by_homology = truekeep_flag = keep
  • 假基因分流:
    • 有明确同源支持但存在 frameshift 或提前终止子,可 is_pseudogene = true,记录到 pseudogene.gff3,不进入高可信蛋白编码基因集合,但可保留在整体 GFF 中作注释。

高可信基因列表 gene.filtered_gene_ids.txt 为后续 07–10 的唯一蛋白编码基因入口。


5.9 07_functional_homology.py —— 同源注释

目的
利用 DIAMOND 对高可信蛋白进行同源比对,给出精简 best hit 表。

输入

  • pep.filtered.fa(06 步从高可信基因提取)
  • db.swissprot、可选 db.nr

输出

  • results/07_homology/homology.best.tsv

表头约定:

1
2
3
4
5
6
7
8
9
10
11
gene_id
protein_id # 如有独立蛋白 ID
db_source # swissprot / nr / other
best_hit_id
best_hit_desc
identity
query_coverage
target_coverage
evalue
bit_score
species # 命中物种简要

5.10 08_functional_interpro.py —— 结构域与 GO 注释

目的
通过 InterProScan 探测结构域与 GO/Pathway 信息。

输入

  • pep.filtered.fa
  • db.interproscan

输出

  • results/08_interpro/interpro.domains.tsv
  • 可选:results/08_interpro/interpro.by_gene.tsv

表头约定:

interpro.domains.tsv

1
2
3
4
5
6
7
8
9
10
gene_id
protein_id
domain_id # Pfam/SMART 等
domain_name
start_aa
end_aa
evalue
go_terms # 以 ; 分隔
pathways # 以 ; 分隔
source_db # pfam / smart / prosite 等

interpro.by_gene.tsv(可选汇总):

1
2
3
4
5
gene_id
domains # "PF00001|PF00002|..."
go_terms
pathways
n_domains

5.11 09_functional_kegg.py —— KEGG 注释(本地 or 线上)

目的
为高可信基因赋予 KO 与 Pathway 信息。

输入:

  • 本地模式:
    • pep.filtered.fa + db.kegg
  • 线上模式:
    • functional.external_kegg_tsv(KAAS 输出)

统一输出:

  • results/09_kegg/kegg.tsv

表头约定:

1
2
3
4
5
6
gene_id
ko # KO 号,如 K00001
kegg_pathways # 以 ; 分隔的 pathway ID,如 ko00010;ko00020
kegg_modules # 可选,; 分隔
definition # KEGG KO 描述(可选)
source # local / kaas

5.12 10_functional_eggnog.py —— eggNOG 注释(本地 or 线上)

目的
通过 eggNOG ortholog group 提供 OG、COG、GO、KO 等。

输入

  • 本地 eggNOG DB 或 functional.external_eggnog_tsv

输出

  • results/10_eggnog/eggnog.tsv

表头约定:

1
2
3
4
5
6
7
gene_id
og # ortholog group ID
cog_category # 如 J, K, L 等
go_terms # 以 ; 分隔
ko # 以 ; 分隔
egg_name # eggNOG 功能名
egg_desc # 功能描述

5.13 10b_ncrna_annot.py —— ncRNA 注释(可选)

目的
注释 tRNA、rRNA、snRNA、snoRNA、miRNA 等 ncRNA。

输入

  • paths.genome.clean 或 masked genome
  • db.rfam
  • tRNAscan-SE 二进制

输出

  • results/10b_ncrna/ncrna.raw.gff3
  • results/10b_ncrna/ncrna.filtered.gff3
  • results/10b_ncrna/TS.ncrna.gff3
  • results/10b_ncrna/ncrna.stats.tsv

表头约定:

ncrna.stats.tsv

1
2
3
4
5
6
rna_class        # tRNA / rRNA / snRNA / snoRNA / miRNA ...
n_loci
median_length_nt
length_min_nt
length_max_nt
note

TS.ncrna.gff3attributes 推荐包含:

  • ID
  • Name
  • biotype

5.14 11_build_publish_release.py —— 发布集构建

目的
整合结构注释、功能注释和 QC 结果,构建统一发布集。

结构输入

  • 高可信蛋白编码基因 GFF:genes.filtered.gff3
  • cds.filtered.fa
  • pep.filtered.fa
  • 可选 ncRNA GFF:TS.ncrna.gff3

功能输入

  • homology.best.tsv(07)
  • interpro.by_gene.tsv(08)
  • kegg.tsv(09)
  • eggnog.tsv(10)
  • gene.confidence.tsv(06)——用于标记 TE-like / pseudogene / novel 等

QC 输入

  • genome BUSCO 摘要
  • proteome BUSCO 摘要
  • gene.stats.tsvaed_distribution.tsv

输出

  • TS.genome.fa
  • TS.genes.gff3
  • TS.cds.fa
  • TS.pep.fa
  • TS.tx2gene.tsv
  • TS.functional.tsv
  • TS.busco_proteome_summary.txt
  • TS.qc_summary.tsv
  • README.annotation.txt

表头约定:

TS.tx2gene.tsv

1
2
tx_id
gene_id

TS.functional.tsv(统一功能总表,推荐结构):

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
gene_id
gene_symbol # 可选
biotype # protein_coding / pseudogene / lncRNA 等(ncRNA 不一定全部在此表)
primary_desc # 优先来源于 SwissProt 或 eggNOG
primary_source # swissprot / eggnog / other
ko # 以 ; 分隔(来自 KEGG + eggNOG)
kegg_pathways # 以 ; 分隔
go_bp # biological process,; 分隔
go_cc # cellular component
go_mf # molecular function
interpro_domains # "PF00001|PF00002|..."
egg_og # eggNOG ortholog group
egg_cog_category # eggNOG/COG 类别
is_te_like # true/false
is_pseudogene # true/false
is_novel_candidate # true/false(12 步)
confidence_level # high/medium/low(来自 gene.confidence)
notes # 其他补充

TS.qc_summary.tsv

1
2
3
metric
value
note

5.15 12_novel_gene.py —— 新基因候选(可选)

目的
在功能注释基础上,筛选缺乏同源或结构域但具表达证据、AED 良好的潜在新基因。

输入

  • TS.functional.tsv
  • gene.confidence.tsv

输出

  • results/12_novel_gene/novel_candidates.tsv

表头约定:

1
2
3
4
5
6
7
8
9
gene_id
has_homology # true/false
has_domain # true/false
expression_support # 如 high/medium/low
aed
is_te_like
is_pseudogene
novel_score # 综合评分
keep_as_novel # true/false

5.16 12b_comparative_qc.py —— 交叉流水线与近缘物种对照 QC(可选)

目的
从两个方向对注释结果进行高级 QC:

  1. 与其它注释流水线 / 外部 GFF 对照;
  2. 与近缘物种注释特征对照。

输入

  1. 本工程注释:
    • TS.genes.gff3
    • TS.pep.fa
  2. 其它流水线注释(可选):
    • comparative_qc.extra_gffs 列出的 GFF
  3. 近缘物种注释:
    • comparative_qc.related_species[*].gff
    • comparative_qc.related_species[*].pep

输出

  • results/12b_comparative_qc/gene_feature_comparison.tsv
  • results/12b_comparative_qc/pipeline_diff_summary.tsv
  • 可选:可视化辅助数据表

表头约定:

gene_feature_comparison.tsv(跨物种结构特征比较):

1
2
3
4
5
6
7
8
9
10
11
12
13
14
species
source # gene-anno / public / other_pipeline
n_genes
mean_cds_len_nt
median_cds_len_nt
mean_exon_num
median_exon_num
single_exon_ratio
mean_intron_len_nt
median_intron_len_nt
busco_c
busco_f
busco_m
note

pipeline_diff_summary.tsv(多流水线对同一物种的注释差异):

1
2
3
4
5
6
7
locus_id             # 标准化 locus 标识(可内部统一)
pipeline_1_gene_id
pipeline_2_gene_id
comparison_type # identical / structure_diff / only_pipeline_1 / only_pipeline_2
delta_exon_num
delta_cds_len_nt
note

6. 功能注释整合机制(本地 vs 线上)

6.1 模式控制

通过 functional.use_local_kegg_dbfunctional.use_local_eggnog_db 控制:

  • true:跑本地数据库;
  • false:使用 functional.external_*_tsv

6.2 合并原则(蛋白编码基因)

在 11 步中:

  1. gene_id 为主键左连接:
    • homology.best.tsv
    • interpro.by_gene.tsv
    • kegg.tsv
    • eggnog.tsv
    • gene.confidence.tsv
  2. primary_desc
    • 优先采用 SwissProt 描述;
    • 无 SwissProt 时,采用 eggNOG/NR 描述。
  3. ko
    • 合并 KEGG + eggNOG 的 KO,去重。
  4. go_bp/cc/mf
    • 合并 InterPro 和 eggNOG 的 GO,拆分到 BP/CC/MF 三列。
  5. interpro_domains
    • 由 InterPro 汇总。
  6. is_te_like / is_pseudogene / is_novel_candidate / confidence_level
    • 直接从 gene.confidence.tsvnovel_candidates.tsv 映射。

6.3 ncRNA 的整合方式

  • ncRNA 的功能信息相对简单,主要通过 biotypeName 表达;
  • 推荐使用独立 GFF:TS.ncrna.gff3
  • 是否合并入 TS.genes.gff3 由项目决定;
  • TS.functional.tsv 不强制包含 ncRNA 信息,可后续视需要扩展 ncrna.functional.tsv

6.4 对下游的意义

  • phylo / aphylo / rna-seq 只需识别:
    • 蛋白编码基因功能表:TS.functional.tsv
    • 结构注释:TS.genes.gff3
    • ncRNA:TS.ncrna.gff3(如需要)。
  • 下游无需关心 KEGG / eggNOG 的具体来源(本地 / 线上)。

7. 日志规范(文件 + 终端流式输出)

7.1 日志文件命名

  • 每个步骤脚本生成一个主日志文件:
    • logs/01_genome_qc.log
    • logs/02_repeat_annot.log
    • ……
  • 如步骤内部有子任务,可在日志中以清晰前缀区分,而不是拆分多个日志文件。

7.2 日志内容要求

每个日志至少包含:

  1. 运行基本信息
    • 起止时间
    • 所用脚本名与简要版本(如内部版本号)
    • 环境信息(可选:hostname、线程数等)
  2. 配置摘要
    • 关键路径(输入 genome、输出目录等)
    • 用到的数据库路径
    • 主参数(如最小 CDS 长度、TE 阈值、是否启用长读模块等)
  3. 输入摘要
    • 样本数 / 记录数 / 基因数等
    • 文件是否存在的检查结果
  4. 核心统计
    • 每一步的关键 QC 数字(如 BUSCO 指标、mapping rate、n_genes 等)
    • 06 步建议记录高可信基因数、假基因数、TE-like 基因数、AED<0.5 的比例等。
  5. 告警与错误
    • WARNING:参数略偏极端、部分输入为空、某些基因被大规模过滤等;
    • ERROR:外部程序返回非零状态、必需文件缺失、表头不匹配等。

7.3 终端流式输出

  • 所有关键信息既写入日志文件,也通过 logger 的 stream 处理器输出到终端;
  • 运行时默认可在终端实时观察进度和关键统计;
  • 有条件时可加上简单的进度计数(如“处理到第 x / N 个样本”)。

8. QC 与交付物

8.1 按阶段的 QC 要点

  1. 基因组层面(01)
    • BUSCO(genome) 完整度
    • N50、GC%、contig 数
  2. 重复注释(02)
    • TE 占比
    • masked genome 与原始 genome 长度差
  3. RNA-seq(03–04)
    • clean reads 数量与质量分布
    • mapping rate / uniquely mapping rate
  4. 长读支持(03b,可选)
    • 长读转录本长度和外显子数分布
    • 支持的剪接位点数量
    • 与短读 hints 的重合度
  5. 基因集 QC(05–06)
    • 基因总数、平均 CDS 长度、平均外显子数
    • 单外显子基因比例
    • 与 TE 重叠的基因比例
    • BUSCO(proteome) 完整度
    • AED 分布(如 >90% 基因 AED < 0.5)
    • 高可信基因数、假基因数、TE-like 基因数等
    • 白名单回捞后被保留的短基因数量
  6. 功能注释(07–10)
    • 有同源注释的基因比例
    • 有 InterPro 域的基因比例
    • 有 KEGG KO 的基因比例
    • 有 eggNOG 注释的基因比例
  7. ncRNA 注释(10b,可选)
    • 各类 ncRNA 的数量及长度分布
    • 与近缘物种的 tRNA/rRNA 数量对比
  8. 发布集与高级 QC(11–12b)
    • 发布集蛋白 BUSCO 指标
    • GFF/CDS/pep/tx2gene 一致性(可在 11 步检查)
    • 与近缘物种的基因结构特征对比(12b 输出)
    • 与其他流水线结果的差异定位(12b 输出)

8.2 对外发布物清单

results/11_publish/ 至少包含:

  • TS.genome.fa
  • TS.genes.gff3
  • TS.cds.fa
  • TS.pep.fa
  • TS.tx2gene.tsv
  • TS.functional.tsv
  • TS.busco_proteome_summary.txt
  • TS.qc_summary.tsv
  • README.annotation.txt

如启用扩展:

  • TS.ncrna.gff3
  • results/12b_comparative_qc/*
  • results/06_gene_qc/aed_distribution.tsv
  • results/06_gene_qc/gene.confidence.tsv(内部 QC 也可在必要时作为附件提供)

9. 方法学模板(占位说明)

本章只保留结构框架,不在蓝图中放完整中英文 Methods 文本。

  • 9.1 中文方法学结构建议
    • 基因组数据来源与质量评估
    • 重复序列注释
    • RNA-seq 数据质控与比对
    • 长读转录本支持(如适用)
    • 基因结构预测与证据整合
    • 基因集质量控制与高可信基因筛选(含 TE 过滤、证据打分与 AED)
    • 功能注释:同源、结构域、KEGG、eggNOG
    • 非编码 RNA 注释(如适用)
    • 注释发布集构建与 QC 汇总
    • 比较质量评估与跨物种对照(如适用)
  • 9.2 英文 Methods 结构建议
    • Genome data and quality assessment
    • Repeat annotation and genome masking
    • RNA-seq processing and genome alignment
    • Long-read transcriptome support (if applicable)
    • Gene prediction and evidence integration
    • Gene set quality control and high-confidence selection (including TE filtering, evidence scoring and AED)
    • Functional annotation (homology, domains, KEGG, eggNOG)
    • Non-coding RNA annotation (if applicable)
    • Construction of the genome annotation release
    • Comparative quality assessment and cross-species comparison (if applicable)

完整中英方法学文本在论文阶段单独撰写,不附于蓝图。


10. 施工注意事项(硬件与数据组织)

  1. 路径规范
    • 一律使用相对路径,通过 config.yaml 管理;
    • 脚本中不写死绝对路径。
  2. 数据库与磁盘
    • BUSCO、InterPro、NR、eggNOG、KEGG 等数据库建议放在 SSD/NVMe 上;
    • 多项目共享同一 db/ 目录;
    • 长读数据体积大时注意 I/O 与临时目录容量。
  3. 脚本与蓝图一致性
    • 调整脚本名称或模块时,要同步修改蓝图中对应章节;
    • 蓝图中提到的每个脚本、目录、表头,都应在工程中可查。
  4. 多服务器场景
    • 核心数据库尽量本地部署;
    • KEGG / eggNOG 部署成本高的服务器可以只读线上结果表(KAAS / eggNOG-mapper);
    • 长读分析可以在专门机器完成,再把 BAM/GTF 作为 03b 输入。
  5. 表头与契约管理
    • 本蓝图中给出的所有 TSV 表头应视为“契约接口”;
    • 修改表头时需同步更新蓝图与下游脚本;
    • 表头变化属于工程级变动,应在 README 中记录。
  6. 日志与 QC 文档化
    • 关键 QC 指标建议整理进 TS.qc_summary.tsv 和 README 中;
    • 对于重要物种,可将 gene_feature_comparison.tsvaed_distribution.tsv 等作为论文补充材料备选。

附录

config.yaml

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
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
# =====================================================================
# gene-anno / config/config.yaml
# 适配《gene-anno 工程蓝图》的完整配置文件示例(以大竹蛏为例)
# 说明:
# - 所有路径推荐使用相对路径
# - 物种与项目信息集中在 project / species 区块
# - 资源(线程、内存)集中在 resources 区块
# - 各步骤(01–12 + 可选模块)通过对应 section 控制
# =====================================================================

# ---------------------------------------------------------------------
# 1. 项目与物种信息(与发布集命名、日志、下游工程关联)
# ---------------------------------------------------------------------
project:
name: "gene-anno" # 工程名称,用于日志、标签等
genome_assembly_name: "DZC_assembly_v1" # 组装版本名称(可用于 Methods/README)
annotation_release_name: "DZC_gene_anno_v1" # 注释版本名称

species:
code: "DZC" # 物种代码(会体现在文件前缀中)
name: "Sinonovacula constricta" # 学名
common_name: "大竹蛏" # 中文名/俗名
taxon_id: 210451 # NCBI Taxon ID(无则可暂写 0)

# ---------------------------------------------------------------------
# 2. 路径配置(全部相对路径,数据与主要结果位置)
# ---------------------------------------------------------------------
paths:
# 基因组序列
genome:
raw: "data/genome/DZC.genome.raw.fa" # 原始基因组(可选备份)
clean: "data/genome/DZC.genome.clean.fa" # 清洗后的基因组(主输入)

# 重复序列相关(02_repeat)
repeat:
masked_genome: "results/02_repeat/DZC.genome.masked.fa" # 屏蔽重复后的基因组
library: "results/02_repeat/DZC.repeat.lib.fa" # 物种特异重复库
gff: "results/02_repeat/DZC.repeats.gff3" # 重复注释 GFF3

# RNA-seq 原始 / clean / 比对结果(03 + 04)
rnaseq:
meta: "data/rnaseq/samples.meta.tsv" # RNA-seq 样本信息表
raw_dir: "data/rnaseq/raw" # 原始 fq.gz
clean_dir: "data/rnaseq/clean" # 质控后的 clean fq
qc_summary: "results/03_rnaseq_qc/fastp_summary.tsv" # fastp 汇总表
mapping_dir: "results/04_mapping" # 比对结果根目录
mapping_bam: "results/04_mapping/rnaseq.sorted.bam" # 合并/汇总后的 BAM(如采用)
mapping_hints_gff: "results/04_mapping/hints.gff" # 由短读生成的 splice hints
mapping_stats_tsv: "results/04_mapping/mapping.stats.tsv"

# 长读 / Iso-Seq / Nanopore 相关路径(03b)
longread:
raw_dir: "data/longread/raw" # 原始长读数据存放位置(可选)
bam: "data/longread/DZC.longread.genome.bam" # 长读对基因组比对后的 BAM
gtf: "data/longread/DZC.longread.transcripts.gtf" # 长读转录本 GTF
hints_gff: "results/03b_longread_support/longread_hints.gff"
stats_tsv: "results/03b_longread_support/longread_transcripts.stats.tsv"

# 近缘物种蛋白(05 步作为外源蛋白证据)
proteins:
dir: "data/proteins" # 近缘物种蛋白所在目录
merged: "data/proteins/homologs.merged.pep.fa" # 合并后的参考蛋白集

# BUSCO 输出目录(01 & 06)
qc:
genome_busco_dir: "results/01_genome_qc/genome.busco" # 基因组 BUSCO 输出目录
genome_stats_tsv: "results/01_genome_qc/genome.stats.tsv"
proteome_busco_dir: "results/06_gene_qc/proteome.busco" # 蛋白 BUSCO 输出目录

# 功能注释结果根目录(07–10)
functional:
homology_dir: "results/07_homology" # 同源注释结果目录
interpro_dir: "results/08_interpro" # InterProScan 结果目录
kegg_dir: "results/09_kegg" # KEGG 注释结果目录
eggnog_dir: "results/10_eggnog" # eggNOG 注释结果目录

# ncRNA 注释输出路径(10b,可选)
ncrna:
out_dir: "results/10b_ncrna" # ncRNA 注释结果目录
raw_gff3: "results/10b_ncrna/ncrna.raw.gff3" # 合并的原始 ncRNA 命中
filt_gff3: "results/10b_ncrna/ncrna.filtered.gff3" # 过滤后的 ncRNA
final_gff3: "results/10b_ncrna/DZC.ncrna.gff3" # 最终发布用 ncRNA 注释
stats_tsv: "results/10b_ncrna/ncrna.stats.tsv"

# 发布集(11_publish)
publish:
out_dir: "results/11_publish" # 发布集总目录
genome_fa: "results/11_publish/DZC.genome.fa" # 发布用基因组
genes_gff3: "results/11_publish/DZC.genes.gff3" # 基因结构注释 GFF
cds_fa: "results/11_publish/DZC.cds.fa" # CDS 序列
pep_fa: "results/11_publish/DZC.pep.fa" # 蛋白序列
tx2gene_tsv: "results/11_publish/DZC.tx2gene.tsv" # 转录本-基因映射
functional_tsv: "results/11_publish/DZC.functional.tsv" # 功能注释总表
proteome_busco_summary: "results/11_publish/DZC.busco_proteome_summary.txt"
qc_summary_tsv: "results/11_publish/DZC.qc_summary.tsv" # QC 汇总表
readme: "results/11_publish/README.annotation.txt" # 发布集说明

# 高级比较 QC 模块 12b 输出路径
comparative_qc:
out_dir: "results/12b_comparative_qc" # 对照 QC 结果目录
gene_feature_comparison_tsv: "results/12b_comparative_qc/gene_feature_comparison.tsv"
pipeline_diff_summary_tsv: "results/12b_comparative_qc/pipeline_diff_summary.tsv"

# ---------------------------------------------------------------------
# 3. 数据库路径(与蓝图 4.2 对应)
# ---------------------------------------------------------------------
db:
busco_lineage: "db/busco/mollusca_odb10" # BUSCO 谱系库
dfam: "db/dfam" # TE 库(Dfam 或 EDTA 结果)
rfam: "db/rfam" # Rfam CM 库,用于 ncRNA
interproscan: "db/interpro" # InterProScan 数据库
nr: "db/nr" # NCBI NR(全库)
swissprot: "db/swissprot" # SwissProt(精库)
kegg: "db/kegg" # 本地 KEGG 数据库(如部署)
eggnog: "db/eggnog" # 本地 eggNOG 数据库

# ---------------------------------------------------------------------
# 4. 功能注释相关设置(07–10,与蓝图第 6 章对应)
# ---------------------------------------------------------------------
functional:
# KEGG / eggNOG 模式切换:
# - true:使用本地数据库并在本地运行
# - false:不跑本地,读取 external_*_tsv
use_local_kegg_db: true
use_local_eggnog_db: true

# 线上 KAAS / eggNOG-mapper 输出文件(当 use_local_* 为 false 时生效)
external_kegg_tsv: "data/functional_external/kegg_online.tsv"
external_eggnog_tsv: "data/functional_external/eggnog_online.tsv"

# 07_functional_homology 参数与输出
homology:
db_priority: ["swissprot", "nr"] # primary_desc 的优先级顺序
max_target_seqs: 5 # 每条蛋白最多保留多少 hits
evalue: "1e-5"
min_identity: 0.30 # 最低序列一致性(0-1)
min_query_coverage: 0.50 # 最低 query 覆盖度
min_target_coverage: 0.50 # 最低 target 覆盖度
best_tsv: "results/07_homology/homology.best.tsv"

# 08_functional_interpro 参数与输出
interpro:
# InterProScan 运行模式配置(可按需调整)
applications: ["Pfam", "SMART", "PROSITEProfiles", "PANTHER"]
use_goterms: true
use_pathways: true
domtblout_tsv: "results/08_interpro/interpro.domains.tsv"
by_gene_tsv: "results/08_interpro/interpro.by_gene.tsv"

# 09_functional_kegg 参数与输出
kegg:
max_hits_per_gene: 5
evalue: "1e-5"
# 输出标准化表(蓝图中 kegg.tsv 的契约)
output_tsv: "results/09_kegg/kegg.tsv"

# 10_functional_eggnog 参数与输出
eggnog:
min_og_support: 1 # 至少多少 ortholog 支持
allow_multiple_ko: true # 是否允许同一基因多个 KO
# 输出标准化表(蓝图中 eggnog.tsv 的契约)
output_tsv: "results/10_eggnog/eggnog.tsv"

# 11 步功能整合相关的一些策略(蓝图 6.2)
merge_policy:
primary_desc_priority: ["swissprot", "eggnog", "nr"] # primary_desc 来源优先级
merge_ko_from: ["kegg", "eggnog"] # ko 列从哪些来源合并
deduplicate_go: true # 合并不同来源 GO 时是否去重
split_go_by_ontology: true # 拆分到 go_bp/go_cc/go_mf 三列

# ---------------------------------------------------------------------
# 5. RNA-seq 样本信息(meta 表)约束(蓝图 5.4)
# ---------------------------------------------------------------------
rnaseq_meta:
file: "data/rnaseq/samples.meta.tsv" # 样本信息表路径
required_columns: # 运行前会检查这些列是否存在
- sample_id
- run_accession
- library_layout
- read_length
- tissue
- condition
- replicate

# ---------------------------------------------------------------------
# 6. BUSCO 配置(01 genome + 06 proteome)
# ---------------------------------------------------------------------
busco:
program: "busco" # BUSCO 可执行文件名
lineage: "db/busco/mollusca_odb10" # 谱系库路径
mode:
genome: "genome" # 基因组模式
proteome: "proteins" # 蛋白模式
threads: 16

# ---------------------------------------------------------------------
# 7. BRAKER3 配置(05_braker3_run)
# ---------------------------------------------------------------------
braker:
program: "braker.pl" # BRAKER 主程序
use_rnaseq: true # 是否使用 RNA-seq 证据
use_proteins: true # 是否使用蛋白证据

augustus_species_name: "generic" # AUGUSTUS 物种模型(可替换为训练后的专用名称)
working_dir: "results/05_braker" # BRAKER 运行目录

# 短读 RNA-seq 证据(来自 04)
rnaseq_bam: "results/04_mapping/rnaseq.sorted.bam"
hints_gff: "results/04_mapping/hints.gff"

# 长读 hints(来自 03b,可选)
longread_hints_gff: "results/03b_longread_support/longread_hints.gff"

# 用作额外证据的近缘物种蛋白
protein_files:
- "data/proteins/homologs.merged.pep.fa"

threads: 32
extra_args: "" # 需要额外传给 BRAKER 的参数可写在这里

# BRaker 主要输出文件(方便后续脚本统一引用)
outputs:
gff3: "results/05_braker/braker.gff3"
cds_fa: "results/05_braker/braker.cds.fa"
pep_fa: "results/05_braker/braker.proteins.fa"
log_summary_tsv: "results/05_braker/braker.log.summary.tsv"

# ---------------------------------------------------------------------
# 8. 基因集 QC 规则(06_gene_qc_and_filter,与蓝图 5.8 完全契合)
# ---------------------------------------------------------------------
gene_qc:
out_dir: "results/06_gene_qc"
gene_stats_tsv: "results/06_gene_qc/gene.stats.tsv"
gene_confidence_tsv: "results/06_gene_qc/gene.confidence.tsv"
high_conf_gene_list: "results/06_gene_qc/gene.filtered_gene_ids.txt"
pseudogene_gff3: "results/06_gene_qc/DZC.pseudogene.gff3"
aed_distribution_tsv: "results/06_gene_qc/aed_distribution.tsv"
proteome_busco_dir: "results/06_gene_qc/proteome.busco"

# 长度与 N 碱基过滤
min_cds_len_nt: 150 # 最小 CDS 长度(碱基数)
min_protein_len_aa: 50 # 最小蛋白长度(氨基酸数)
max_ambiguous_nt: 0 # 允许的 N 碱基数量(>0 时视为低可信)

# TE 相关过滤
filter_te_overlap: true
te_overlap_gff: "results/02_repeat/DZC.repeats.gff3"
te_overlap_fraction: 0.5 # CDS 与 TE 重叠比例超过该值则标记为 TE-like

# AED 相关
aed_source:
use_rnaseq: true # 是否利用 RNA-seq 计算 AED
use_longread: true # 是否利用长读计算 AED(enable=true 时)
aed_good_threshold: 0.5 # AED <= 该值视为拟合良好
aed_bins: # 生成 aed_distribution.tsv 时使用的分箱
- [0.0, 0.1]
- [0.1, 0.2]
- [0.2, 0.3]
- [0.3, 0.4]
- [0.4, 0.5]
- [0.5, 0.6]
- [0.6, 0.7]
- [0.7, 0.8]
- [0.8, 0.9]
- [0.9, 1.0]

# BUSCO 相关
keep_partial_busco: true # 部分 BUSCO 是否保留
busco_summary_file: "results/06_gene_qc/proteome.busco/short_summary.txt"

# 证据打分(E-score)配置 —— 蓝图中的 evidence_score 来源与权重示例
evidence_scoring:
# 每项证据的加分(可按项目需要调整)
weights:
rnaseq_support: 2 # has_rnaseq_support = true
longread_support: 2 # has_longread_support = true
swissprot_hit: 3 # has_swissprot_hit = true
interpro_domain: 2 # has_interpro_domain = true
eggnog_og: 1 # has_eggnog_og = true
busco_complete: 3 # busco_status = C
busco_fragmented: 1 # busco_status = F
low_aed: 2 # aed <= aed_good_threshold
high_te_penalty: -3 # is_te_like = true
short_cds_penalty: -2 # is_short_cds = true
frameshift_penalty: -3 # has_frame_shift_or_stop = true

# 基于 E-score 的置信度分级阈值(仅示例,可按实际调整)
confidence_thresholds:
high: 5 # evidence_score >= 5 → high
low: 1 # evidence_score <= 1 → low,其余为 medium

# 白名单回捞与假基因分流策略
rescue_policy:
# 满足任意一条即视为“白名单,可回捞”
require_swissprot_or_domain: true # SW 或 InterPro 任一存在即可回捞
min_rescue_score: 2 # 即使结构较差,只要 E-score >= 2 即可考虑回捞

pseudogene_policy:
# 判定假基因的条件(需配合脚本中逻辑使用)
require_homology: true # 仅对有同源注释的基因评估假基因
frameshift_required: true # 有 frameshift 或提前终止子
allow_keep_in_gff: true # 是否在整体 GFF 中保留假基因记录

# ---------------------------------------------------------------------
# 9. 长读转录本支持模块配置(03b_longread_support,可选)
# ---------------------------------------------------------------------
longread_support:
enable: false # 是否启用长读支持模块(03b),无长读时为 false
bam: "data/longread/DZC.longread.genome.bam" # 长读比对到基因组后的 BAM
gtf: "data/longread/DZC.longread.transcripts.gtf" # 长读转录本 GTF(如 IsoSeq3/TAMA 结果)
assembler: "" # 长读转录本流程名称(如 "isoseq3", "tama", "stringtie2")
min_supporting_reads: 2 # 生成 hints 时,最少支持 reads 数
min_isoform_length_nt: 300 # 过滤过短转录本的阈值(碱基)
hints_gff: "results/03b_longread_support/longread_hints.gff"
stats_tsv: "results/03b_longread_support/longread_transcripts.stats.tsv"

# ---------------------------------------------------------------------
# 10. ncRNA 注释模块配置(10b_ncrna_annot,可选)
# ---------------------------------------------------------------------
ncrna:
enable: false # 是否启用 ncRNA 注释模块(10b)
rfam_cm: "db/rfam/Rfam.cm" # Rfam CM 数据库文件(需事先 cmpress)
cmscan_program: "cmscan" # Infernal cmscan 可执行文件
trnascan_program: "tRNAscan-SE" # tRNAscan-SE 可执行文件
score_cutoff: 20.0 # cmscan 命中得分的最低阈值(示例)
report_classes: # 希望在最终 GFF 中保留的 ncRNA 类型
- "tRNA"
- "rRNA"
- "snRNA"
- "snoRNA"
- "miRNA"
outputs:
raw_gff3: "results/10b_ncrna/ncrna.raw.gff3"
filt_gff3: "results/10b_ncrna/ncrna.filtered.gff3"
final_gff3: "results/10b_ncrna/DZC.ncrna.gff3"
stats_tsv: "results/10b_ncrna/ncrna.stats.tsv"

# ---------------------------------------------------------------------
# 11. 新基因候选模块配置(12_novel_gene,可选)
# ---------------------------------------------------------------------
novel_gene:
enable: false # 是否启用新基因候选筛选(12)
input_functional_tsv: "results/11_publish/DZC.functional.tsv"
input_gene_confidence_tsv: "results/06_gene_qc/gene.confidence.tsv"
novel_candidates_tsv: "results/12_novel_gene/novel_candidates.tsv"

# 新基因候选判定逻辑(结合蓝图 5.15)
criteria:
require_no_homology: true # has_homology = false
require_no_domain: true # has_domain = false
max_aed: 0.4 # AED <= 该值
disallow_te_like: true # is_te_like = false
disallow_pseudogene: true # is_pseudogene = false
min_expression_support: "medium" # expression_support 至少达到 medium

scoring:
base_score: 0
add_expression_high: 3
add_expression_medium: 1
penalty_te_like: -5
penalty_low_confidence: -2
keep_threshold: 2 # novel_score >= 2 → keep_as_novel = true

# ---------------------------------------------------------------------
# 12. 高级比较 QC 模块配置(12b_comparative_qc,可选)
# ---------------------------------------------------------------------
comparative_qc:
enable: false # 是否启用 12b 比较 QC 模块

# 其它注释流水线的 GFF,用于与 BRAKER3 结果对照(可为空列表)
extra_gffs: [] # 例如: ["data/annotation/other_pipeline1.gff3"]

# 近缘物种列表,用于基因集结构特征对比
related_species:
- name: "RefSpecies1" # 近缘物种名称
gff: "data/annotation/ref1.genes.gff3" # 近缘物种 GFF
pep: "data/annotation/ref1.pep.fa" # 近缘物种蛋白序列
# - name: "RefSpecies2"
# gff: "data/annotation/ref2.genes.gff3"
# pep: "data/annotation/ref2.pep.fa"

outputs:
gene_feature_comparison_tsv: "results/12b_comparative_qc/gene_feature_comparison.tsv"
pipeline_diff_summary_tsv: "results/12b_comparative_qc/pipeline_diff_summary.tsv"

# ---------------------------------------------------------------------
# 13. 资源配置(线程与内存)
# ---------------------------------------------------------------------
resources:
threads:
genome_qc: 16
repeat_annot: 16
rnaseq_qc: 8
rnaseq_mapping: 16
longread_support: 8
braker3: 32
gene_qc: 8
functional_homology: 16
functional_interpro: 32
functional_kegg: 8
functional_eggnog: 8
ncrna_annot: 8
publish: 4
novel_gene: 4
comparative_qc: 4

memory_gb:
genome_qc: 32
repeat_annot: 64
rnaseq_qc: 16
rnaseq_mapping: 32
longread_support: 32
braker3: 128
gene_qc: 16
functional_homology: 64
functional_interpro: 128
functional_kegg: 32
functional_eggnog: 32
ncrna_annot: 32
publish: 16
novel_gene: 16
comparative_qc: 16

# ---------------------------------------------------------------------
# 14. 工具可执行文件路径
# ---------------------------------------------------------------------
tools:
python: "python"
fastp: "fastp"
fastqc: "fastqc"
hisat2: "hisat2"
star: "STAR"
samtools: "samtools"
minimap2: "minimap2"

braker3: "braker.pl"
augustus: "augustus"
genemark: "gmes_petap.pl"

diamond: "diamond"
interproscan: "interproscan.sh"
busco: "busco"

repeatmasker: "RepeatMasker"
repeatmodeler: "RepeatModeler"
edta: "EDTA.pl"

cmscan: "cmscan"
trnascanse: "tRNAscan-SE"
bedtools: "bedtools"
gffread: "gffread"

# ---------------------------------------------------------------------
# 15. 日志配置(文件 + 终端流式输出,与蓝图第 7 章一致)
# ---------------------------------------------------------------------
logging:
log_dir: "logs" # 日志目录
log_level: "INFO" # 日志等级:DEBUG / INFO / WARNING / ERROR
format: "[%(asctime)s] [%(levelname)s] %(message)s"
datefmt: "%Y-%m-%d %H:%M:%S"


基因组注释计划蓝图
https://oldstory.cn/2025/12/20/ji_yin_zu_zhu_shi_ji_hua_lan_tu/
作者
Ricardo
发布于
2025年12月20日
许可协议