从NCBI下载到生成进化树用ProkkaRoary完成细菌泛基因组分析的完整实战记录附批量脚本在微生物基因组学研究中泛基因组分析已成为揭示物种内基因多样性的重要手段。本文将手把手带您完成从NCBI数据下载到最终生成系统发育树的全流程特别适合需要在短时间内掌握完整分析链条的研究人员。不同于碎片化的教程我们聚焦于实际科研场景中的连贯操作并提供可直接复用的批量处理脚本。1. 数据获取与预处理1.1 NCBI基因组数据批量下载首先登录NCBI Assembly数据库在搜索框输入目标菌种名称如Bacillus subtilis。建议设置以下筛选条件Assembly level: Complete Genome或ChromosomeRefSeq category: Reference/Representative genome排除contig/scaffold级别的组装勾选所需基因组后使用Download Assemblies功能选择File type: GenBank (.fna)Include: Genomic sequence (FASTA)实用技巧对于大规模下载推荐使用NCBI的datasets命令行工具datasets download genome accession --inputfile accessions.txt --filename ncbi_dataset.zip其中accessions.txt包含需要下载的基因组Accession编号列表。1.2 数据质量控制解压下载的压缩包后建议进行基本质量检查使用seqkit stats快速统计序列信息seqkit stats *.fna genome_stats.txt重点关注序列数量应为1条/基因组平均序列长度应与物种预期基因组大小匹配典型问题处理若发现多个contig的基因组建议重新筛选更完整的组装版本或使用--compliant参数运行Prokka2. 批量基因组注释流程2.1 Prokka环境配置推荐使用conda管理分析环境conda create -n prokka -c bioconda prokka1.14.6 conda activate prokka prokka --setupdb # 关键步骤初始化数据库2.2 自动化注释脚本创建run_prokka.sh脚本实现批量处理#!/bin/bash mkdir -p gff_output for fna in *.fna; do sample$(basename $fna .fna) prokka $fna \ --prefix $sample \ --outdir ${sample}_annot \ --cpus 8 \ --force # 覆盖已有结果 cp ${sample}_annot/*.gff gff_output/ done关键参数说明--cpus根据服务器资源调整--force确保重复运行时覆盖旧结果--addgenes添加gene特征Roary分析必需3. 泛基因组分析核心步骤3.1 Roary安装与配置虽然可通过conda安装但推荐系统级安装确保稳定性sudo apt-get install -y roary mafft fasttree3.2 运行泛基因组分析准备roary_analysis.sh脚本#!/bin/bash roary -f roary_results \ -p 16 \ -e \ # 使用MAFFT进行多序列比对 -v \ # 详细日志 -i 90 \ # 蛋白相似度阈值 gff_output/*.gff输出文件解读gene_presence_absence.csv核心/可变基因矩阵summary_statistics.txt泛基因组大小等统计量core_gene_alignment.aln核心基因多序列比对3.3 核心基因系统发育树构建使用FastTree快速生成进化树FastTree -nt -gtr roary_results/core_gene_alignment.aln core_genome_tree.nwk对于更精确的分析可考虑iqtree -s core_gene_alignment.aln -m MFP -bb 1000 -nt AUTO4. 结果可视化与解读4.1 泛基因组特征分析使用roary_plots.py生成基础图表roary_plots.py roary_results/core_gene_alignment.aln roary_results/gene_presence_absence.csv将生成泛基因组累积曲线pangenome_matrix.png核心基因热图clustered_pangenome_matrix.png4.2 进化树美化推荐使用FigTree进行树形调整或R语言ggtree包library(ggtree) tree - read.tree(core_genome_tree.nwk) ggtree(tree) geom_tiplab() xlim(0, 0.2) geom_nodepoint(colorred, alpha0.5)4.3 高级分析拓展对于深入挖掘可尝试Panaroo改进的泛基因组分析工具pyseer基因-表型关联分析Scoary基因与表型特征关联研究5. 常见问题解决方案Q1: Prokka运行时出现数据库错误prokka --listdb # 检查数据库状态 prokka --cleandb # 清理损坏的数据库 prokka --setupdb # 重新初始化Q2: Roary运行内存不足增加-s参数降低blastp严格度使用--split_paralogs处理旁系同源基因Q3: 进化树分支支持率低增加比对基因数量降低-i参数使用IQ-TREE代替FastTree实际项目中处理30个芽孢杆菌基因组时使用32核服务器完整流程约需6-8小时。建议将中间结果特别是GFF文件妥善保存便于后续重新分析。