ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

Kraken2+Bracken宏基因组物种注释实战:从安装建库到结果解读

Kraken2+Bracken宏基因组物种注释实战:从安装建库到结果解读 最近两个礼拜我帮好几个课题组处理了宏基因组测序数据发现一个现象很多人的流程跑到了组装和分箱阶段才想起来问我这些contig到底是哪个物种的。其实在宏基因组分析管线里物种注释这一步应该更早介入直接对原始测序读段做分类学归属才能在后续分析里不被一堆来源不明的contig拖慢节奏。Kraken2和Bracken这套组合就是我现在最常用的快速物种注释方案从安装到出结果一条路走下来简单、快、结果还靠谱。这篇文章把安装、建库、运行到常见坑位的完整过程都写出来适合刚转型做宏基因组分析、或者被BLAST速度折磨过的朋友参考。我尽量把每个环节为什么这么做也说清楚因为生信工具最大的门槛从来不是命令本身而是你不理解参数背后的逻辑出了问题根本不知道从哪里排查。1. 从BLAST到K-merKraken2的出现背景与定位1.1 传统方法慢在哪里Kraken2靠什么提速我在刚接触生信那会儿做物种注释用的是BLAST把一条read拿去和NR数据库比对一条read往往要等好几秒甚至更久一个宏基因组样本几百万条reads算下来几十个小时起步而且结果里还充斥着大量低相似度的命中你需要自己定阈值去筛麻烦得很。后来有了MEGAN这类工具用BLAST结果做LCA最低共同祖先归类准确度上去了但瓶颈还是卡在比对这一步速度根本没有质变。Kraken2的思路是完全另一条路它不比对它查表。把每条read切成固定长度的k-mer然后在构建好的k-mer数据库里精确检索每个k-mer都能映射到一个分类学节点最后根据一条read上所有k-mer对应的分类节点找到它们共同的祖先作为这条read的分类归属。整个过程只有哈希查找和树遍历没有动态规划没有打分矩阵所以速度能比BLAST快几个数量级。这也是它名字里K的由来——k-mer。1.2 Kraken2和Bracken的分工关系光有Kraken2你能知道一条read属于哪个物种但宏基因组分析还关心另一个问题这些物种的相对丰度是多少Kraken2自身的report文件虽然是按分类层级统计的但统计口径是被分类的reads数量这在很多场景下会低估或高估物种的真实丰度原因后面细说。Bracken就是专门来解决丰度估算问题的。它利用Kraken2对每条read的分类结果通过贝叶斯概率重新分配那些被归到高分类层级比如属、科的reads把它们往下推到物种级别从而给出更符合真实情况的丰度估计。所以现在的主流做法是明确的流水线分工Kraken2负责分类Bracken负责算丰度。两个工具配合使用才能在有哪些物种和各占多少比例两个维度上都拿到可靠答案。1.3 这套方案适合什么场景从我自己的经验看Kraken2Bracken最适合下面这类工作宏基因组样本的初步物种构成调查想快速知道样本里有哪几类微生物临床样本的病原微生物筛查需要快速排除或确认特定物种宿主DNA污染评估比如肿瘤组织测序里混了多少细菌在正式组装之前对reads做一次身份筛查辅助后续分析策略。如果你手头是扩增子数据16S/ITS那有更适合的QIIME 2或者mothur不推荐拿Kraken2硬上因为扩增子读段本身信息量有限k-mer策略的优势发挥不出来。2. 安装Kraken2和Bracken两条路线与版本选择2.1 用Conda/Mamba隔离环境安装Kraken2和Bracken的依赖非常少核心就是C编译好的二进制文件加一个数据库目录所以安装难度在生信工具里属于入门级别。但我强烈不建议在base环境里装因为Kraken2的Python辅助脚本依赖的包版本太老容易和后续装的其他工具起冲突。我的习惯是用mamba单独创建一个环境mamba create -n kraken2 -c bioconda kraken2 bracken conda activate kraken2mamba比conda的依赖解析速度快太多了特别是处理bioconda这种包数量庞大的频道用conda经常卡在Solving environment这一步换了mamba基本几十秒搞定。如果机器上还没有mamba先装一个conda install -n base -c conda-forge mamba -y装完之后验证一下是否成功kraken2 --version bracken --version正常情况下会看到类似于Kraken2 version 2.1.3、Bracken 2.9这样的版本号输出。如果只是报command not found多半是环境没激活或者是conda源没配好换国内镜像源再装一次就行。2.2 源码编译安装的注意事项conda在国内源的选择比较多但如果是在一台完全离线的内网机器上部署或者安全策略不允许用conda下载预编译包那就得走源码编译。Kraken2的源码在GitHub上clone下来之后在src目录下make就行git clone https://github.com/DerrickWood/kraken2.git cd kraken2/src make编译本身很顺一般几分钟结束。编译完的kraken2和kraken2-build、kraken2-inspect这些可执行文件都在src目录里需要自己往PATH里加或者手动复制到/usr/local/bin。这里有个小坑源码编译默认不带上BrackenBracken也要单独去GitHub拉一份源码来编译git clone https://github.com/jenniferlu717/Bracken.git cd Bracken/src makeBracken编译后生成的bracken、bracken-build这些脚本在根目录而不是src目录。另一个需要注意的点是Bracken里有个核心脚本est_abundance.py它依赖Python的几个科学计算库源码方式安装的话记得把numpy、scipy、pandas这些先装好。2.3 版本选择建议Kraken2的版本更新不算频繁选bioconda默认的最新稳定版就好。但要注意的是Kraken2数据库的构建格式和软件版本存在一定的兼容性关系如果你用的是新版软件配老数据库比如2020年以前构建的读取时可能报错。我一般会记录每个环境的软件版本和数据库版本做项目归档时一并提交方便后续复现和分析溯源。Bracken的版本选择要特别留意它必须和Kraken2数据库配套使用因为Bracken在重估丰度时需要读Kraken2数据库里的k-mer分布信息。不同版本的Bracken对数据库目录结构的要求基本一致但如果Kraken2数据库太老Bracken新版本可能解析不了里面的taxonomy信息。建议安装时用bioconda同时指定版本安装确保数据库兼容而不是一个最新一个最老。3. 数据库准备官方库、标准库构建与自定义库数据库是Kraken2的灵魂软件本身只是查询工具数据库的质量直接决定分类结果的准确度。这一步最容易出问题也最需要耐心我分几个场景详细说。3.1 官方预构建数据库的选择对大多数宏基因组用户来说最省事的方式是直接下载官方预构建的数据库不需要自己从NCBI一份一份拉参考基因组。Kraken2官方在AWS上维护了多个版本的预构建索引截至我写这篇博文比较常用的几个是数据库名称内容磁盘占用解压后内存建议minikraken2_v2细菌、古菌、病毒精简版约4GB8GBStandard细菌、古菌、病毒、人源约50GB约70GBStandard-16细菌、古菌、病毒约28GB约50GBPlusPFStandard 质粒/原生动物/真菌约75GB约100GBEuPathDB46真核病原体为主约80GB约120GB挑数据库最核心的考量是你要分析的对象里可能有什么。如果做人体肠道宏基因组Standard库就够用了如果做环境样本尤其是土壤或海洋真菌和原生动物不能忽略那就上PlusPF如果只是快速检测一个临床样本里有没有特定病原体minikraken2都能对付。数据库越大分类越全但内存和磁盘开销也越大完全可以根据实际需要来。下载官方库的方式我推荐直接用wget加断点续传wget -c https://genome-idx.s3.amazonaws.com/kraken/k2_pluspf_20230314.tar.gz tar -xzvf k2_pluspf_20230314.tar.gz -C $DBDIR数据库解压后是一个目录里面有hash.k2d、opts.k2d、taxo.k2d三个文件Kraken2运行时会同时读取这三个文件。关于下载国内直连AWS的速度时好时坏遇到慢的时候不要反复中断重试用-c参数续传或者选择非高峰期比如凌晨下载。部分国内镜像站也会同步这些数据库可以找一下你所在课题组常用的镜像源。这属于网络环境差异没有统一的最优解能下下来就是胜利。3.2 自己构建Standard库的完整流程如果网络条件实在不允许下载几十GB的预构建库或者你就是想要一个定制化的库比如只需要某些特定物种那kraken2-build命令就派上用场了。自建库的全流程是下载taxonomy 下载参考序列 构建索引三步。先建目录并下载分类学信息mkdir -p krakendb kraken2-build --download-taxonomy --db krakendb这一步会从NCBI下载完整的taxdump文件包括nodes.dmp、names.dmp、merged.dmp等。国内访问NCBI一般问题不大但偶尔也会断建议同样用screen或nohup挂后台执行。接着选参考库。Kraken2官方库里把微生物参考基因组分成了细菌bacteria、古菌archaea、病毒viral、人类human、真菌fungi、植物plant等若干个子库你可以按需下载kraken2-build --download-library bacteria --db krakendb kraken2-build --download-library archaea --db krakendb kraken2-build --download-library viral --db krakendb kraken2-build --download-library human --db krakendb下载完所有需要的library之后开始构建索引kraken2-build --build --threads 32 --db krakendb构建索引这一步非常吃内存标准库构建过程中的峰值内存可能超过100GB我之前在一台128GB内存的服务器上跑都微微出汗。构建结束后同样会得到hash.k2d、opts.k2d、taxo.k2d三个文件。自建库最大的Value是可以完全控制参考基因组列表。比如我处理过一个深海沉积物样本里面有很多未培养微生物通用库分类效果很差我就把已报道的相关候选物种基因组全部下载下来加进自定义库分类效果立刻有了明显改善。Kraken2支持往库里追加序列kraken2-build --add-to-library custom_genome.fna --db krakendb kraken2-build --build --threads 32 --db krakendb这里要注意一个细节custom_genome.fna的FASTA头格式里序列ID必须能对应到NCBI的taxonomy ID否则会报错或者被跳过。Kraken2官方文档推荐用如下形式kraken:taxid|561|NC_000913.3 Escherichia coli str. K-12 substr. MG1655其中561是大肠杆菌的taxonomy ID。如果你的参考序列文件里没有这个信息需要先去NCBI的Assembly数据库查对应物种的taxid然后批量重写FASTA头。这个步骤我建议写个简单脚本来做手动改上百条记录太容易出错了。3.3 数据库构建中的常见失败模式kraken2-build最容易挂在两个地方。第一个是下载人类基因组时NCBI的FTP偶尔会拒绝并发连接导致一批代表序列下载失败。重新执行下载命令时已下载的文件会被跳过Kraken2构建脚本有断点续传的设计所以重试是安全的。第二个是--build步骤报错信息里出现Not all downloaded files were used或者Taxonomy ID not found。前者多是因为部分参考基因组的assembly_summary.txt没有正确同步后者是因为自定义库里的某些taxid在nodes.dmp中不存在。解决思路都是先去检查对应的taxid是否已经包含在你下载的taxonomy里如果确实没有可以单独下载taxdump再替换。4. 第一次跑分类Kraken2命令参数与输出格式精读4.1 一个标准的分步式命令行假设你已经准备了一个双端测序样本文件名为sample_R1.fastq.gz和sample_R2.fastq.gz数据库目录是krakendb常见的运行命令可以写成这样kraken2 --db krakendb \ --paired \ --threads 16 \ --confidence 0.2 \ --minimum-base-quality 15 \ --classified-out classified_#.fq \ --unclassified-out unclassified_#.fq \ --output kraken2.output \ --report kraken2.report \ sample_R1.fastq.gz sample_R2.fastq.gz这些参数的用途我逐个说一下因为它们直接决定结果质量。--paired告诉Kraken2输入是双端reads它会将read1和read2的k-mer合并判断如果不加这个参数两条reads会被当成两条独立的序列去分类。--confidence是分类置信度阈值取值范围0到1。它的底层逻辑是一条read里有多少比例的k-mer支持最终的分类结论。默认值是0意味着只要某个分类节点的支持大于0就采用产生的结果里会包含大量低置信度的分类。做宏基因组多样性研究时我一般建议设0.2到0.4能够过滤掉很多假阳性分类但如果是做病原检测想把弱信号也保留下来那可以设置为0靠后续人工复核。--minimum-base-quality是碱基质量过滤阈值默认0。我给测序数据设15相当于把质量很差的碱基对应的k-mer直接丢弃能减少测序错误导致的虚假分类。这个参数对低质量数据特别有效处理长度较短的reads时不要设太高超过20可能会丢失过多信息。--classified-out和--unclassified-out分别输出分类成功和未分类的reads文件名里的#会被替换成1和2对应双端的两条reads。这两个文件很值钱下游做组装或者提取目标物种序列完全依赖它们。输出结果有两个核心文件--output是每条read的详细分类结果--report是汇总的分类统计报告。很多人直接看report但忽略output后面我详细解读两者差异。4.2 逐行读懂分类输出文件与报告文件output文件的每一行对应一条read用制表符分隔成5列。我摘几行实际输出来看C V00534.1 561 351 351:14872 349:8987 355:5793 561:2350 U V00536.1 0 91 0:91 C M34882.1 1280 243 0:243 1280:1869第一列是分类状态C代表Classified成功分类U代表Unclassified未分类。第二列是read的序列ID。第三列是taxid。第四列是read长度。第五列是k-mer的详细匹配结果格式是分类节点:该节点打上的k-mer数冒号分隔空格分隔多个节点。这一列是理解Kraken2分类逻辑的关键也能用来排查为什么某条read被归类到某个奇奇怪怪的物种上。report文件的格式相对更结构化一共6列2.30 100000 100000 U 0 unclassified 97.70 4250000 4250000 - 1 root 97.70 4250000 4150000 D 2 Bacteria ...这6列分别是该分类层级的reads百分比相对总reads、该层级的reads数、该层级下的累积reads数包含所有子层级、分类层级代码、taxid、物种学名。层级代码U代表unclassifiedD是domain域P是phylum门C是class纲O是order目F是family科G是genus属S是species种S1是亚种。层级代码里还有一个-代表root和中间汇总节点。读report的时候要注意第五列的百分比是该层级自身直接分类到的reads占全部reads的比例不是该层级子层级的累积占比。累积占比要看第三列的累积reads数。很多人做柱状图时拿错列导致物种比例加起来超过100%一般都是这个原因。4.3 关于分类结果的一个实战案例观察我曾经跑过一份模拟群落数据里面按已知比例混合了5株细菌和1株真菌。用默认参数跑了Kraken2后report文件里真核那株真菌的reads占比跟理论比例差别很大反而多出来不少低丰度物种的信号。把--confidence提高到0.3之后大部分低丰度杂信号消失了真菌的reads占比也回归到了理论区间附近。这让我养成了一个习惯正式分析前先拿模拟群落或者已知样本测一遍把参数调到合理区间后再跑批量样本避免默认参数一把梭。5. Bracken丰度重估从分类结果到真实比例5.1 为什么Kraken2的reads数不能直接当丰度用你在Kraken2 report里看到物种A占了40%的reads如果直接把这个数字当作物种A的相对丰度在宏基因组研究里可能会被审稿人挑战原因主要有三个。首先不同物种的基因组大小不一样。一个4Mb基因组的细菌和一个40Mb基因组的真菌打碎成相同长度的reads前者的reads数天然只有后者的十分之一不到。reads占比反映的是基因组大小x拷贝数x丰度的复合信号不单纯是丰度。其次基因组内部的重复序列和高保守区域会造成k-mer多映射。比如一条16S rRNA基因序列几乎在所有细菌里都高度相似它产生的k-mer会同时命中多个物种Kraken2的LCA算法会把它归到更高级别的分类节点比如细菌域导致物种水平reads被压制。Bracken解决的就是这个问题。它充分利用Kraken2已经分类到中间节点的reads结合数据库里各物种基因组的k-mer分布特征用贝叶斯模型把reads重新分配到最可能的物种层级最终输出经过校正的丰度估计值。5.2 用bracken-build构建丰度估计的数据库文件Bracken能在物种级别做重估但它需要知道每个物种基因组里各k-mer的出现频次。这个信息就要通过bracken-build提前算好并存成数据库。操作非常简单bracken-build -d krakendb -t 32 -k 35 -l S参数含义-d指定Kraken2数据库目录-t线程数-k指定k-mer长度必须和构建Kraken2数据库时设定的k-mer长度一致。官方库默认是35大多数自建库也是35如果不确定可以去库目录里看opts.k2d里的参数或者直接看构建时的日志-l指定要计算哪个分类层级的分布S代表species也可以填Ggenus、Ffamily等。实际研究中大多数关注物种水平所以填S为主。bracken-build运行时间取决于数据库大小PlusPF库跑起来可能需要1到2小时。运行结束后会在数据库目录下生成一个名为kmer_distributions_35_species.txt这样的文件这就是后续Bracken预估丰度要用到的核心文件。你可能觉得这一步麻烦其实它只需要在换数据库时跑一次同一个数据库下所有样本都不用重复计算所以时间成本完全可接受。5.3 Bracken运行参数与输出解读Bracken主程序读入Kraken2的report文件输出修正后的丰度列表bracken -d krakendb \ -i kraken2.report \ -o bracken.species.txt \ -r 150 \ -l S \ -t 16这里-r 150是read长度很关键。Bracken重估丰度时依赖reads在基因组上的覆盖分布不同read长度的k-mer覆盖率不同所以尽量填真实测序读长而不是随便填一个。如果测序读长是双端150bp填150是PE100的测序填100。填错了也不会报错但结果会有系统偏差。运行结束后得到的bracken.species.txt前两行是命令和数据库信息正式数据从第三行开始。字段分别是name taxonomy_id taxonomy_lvl kraken_assigned_reads added_reads new_est_reads fraction_total_reads Escherichia coli 562 S 120000 80000 200000 0.2134kraken_assigned_readsKraken2直接分类到该物种的reads数added_readsBracken从更高分类层级比如属、科重新分配下来的reads数new_est_reads两者之和这就是Bracken估算的该物种最终reads数fraction_total_reads该物种占所有分类reads的比例。我在实际项目里用最后这个fraction_total_reads列做物种相对丰度柱状图、PCA分析和差异物种筛选。Bracken结果还能直接接上LEfSe、STAMP这类差异分析工具从Kraken2 report到最终图表整条链路非常顺滑。有的朋友会纠结一个问题Bracken输出的物种列表里的reads数是reads的真实归属吗坦白讲这是一个统计估计不是生物学事实。Bracken在算法层面假设了基因组在样本中均匀覆盖但真实样本里有GC偏好、拷贝数变异、近缘物种同时存在等复杂情况丰度估计始终是近似值。所以我的原则是做看板、做趋势分析用Bracken结果做关键菌株的绝对定量还得靠qPCR或者覆盖度分析两类方法互相印证。6. 实战中绕不开的坑内存、数据库、结果异常的排查链路工具本身不难装但真到跑大数据的时候问题就一个一个冒出来了。我把自己踩过的、以及帮人排过的几个典型问题整理一下按现象-排查链路-最终解决的顺序写方便你在出问题时照方抓药。6.1 报错Failed to load hash.k2d或者段错误怎么看这类问题高度疑似数据库文件不完整Kraken2读库时如果发现哈希文件损坏直接抛异常或者以段错误Segmentation fault形式退出。排查链路如下先检查数据库目录的完整性三个关键文件hash.k2d、opts.k2d、taxo.k2d是否都存在看文件大小是否和官方文档一致。可以拿解压时的tar包大小比对或者重新解压一遍校验如果下载过程中用wget发生过中断即便日志显示完成也可能存在静默损坏。我遇到过好多次文件大小看起来没问题但内部数据错位解决办法是删掉重下用-c续传节省时间的概率反而低检查当前kraken2版本是否和数据库构建时的版本跨度很大如果软件版本过旧读新库时也容易出现异常行为。6.2 内存不足怎么办--memory-mapping的妙用Kraken2标准库跑起来大约要吃60GB到80GB内存如果服务器只有64GB很容易直接被OOM杀掉。解决思路有三个方向第一换小库。minikraken2或者Standard-16是低内存机器的常备选择代价是分类覆盖面变窄很多属种归不到物种级别。第二用Kraken2自带的--memory-mapping参数。这个参数让程序不把整个数据库一次性加载进物理内存而是采用内存映射memory-mapping的方式读数据库文件。系统按需从磁盘读入数据可以显著降低峰值内存代价是速度会慢一些因为频繁的磁盘I/O比纯内存慢得多。我自己在64GB机器上跑PlusPF库时加上这个参数后确实能跑完推荐遇到内存瓶颈的人先试这个方案。第三限制线程数。Kraken2每个线程都会申请一部分临时缓冲区线程数从32降到16内存占用能降不少速度也不会明显折损尤其是当数据库文件在机械硬盘上时线程太多反而在抢I/O。6.3 分类结果大量unclassified问题出在哪如果样本里超过80%的reads都没被分类而你的样本确实应该包含大量微生物那就需要排查了。我经历过的场景里最常见的原因有两类。一类是样本本身宿主污染极重比如组织样本或者血样人体reads占了绝大多数而你没有选带human库的数据库版本或者用了minikraken2这种不含人源序列的小库结果人源reads全部进入unclassified。解决方法是换用Standard/PlusPF库或者在下游统计时直接用Bowtie2把宿主reads比对后剔除再做Kraken2分析。另一类是低复杂度序列太多比如接头序列、polyA尾巴、重复序列。Kraken2对完全重复的k-mer会在数据库里找不到匹配从而归为unclassified。这个可以在上游用fastp做质量控制和接头去除按标准宏基因组流程走一遍就基本能消除。如果你想更细粒度地定位可以用kraken2-inspect查看数据库覆盖了哪些物种和它们的k-mer数量。如果数据库里对应物种的k-mer覆盖就很少那很可能本来就是数据库不完备优先换库而不是调参。6.4 分类结果里出现不可能存在的物种有一回我跑肠道菌群样本结果里冒出了不少海洋细菌第一反应是数据库污染。排查后发现是这些海洋细菌的k-mer在数据库里和其他常见肠道菌高度重叠低复杂度k-mer导致交叉映射。解决办法是提高--confidence阈值并且跑完Bracken后再过滤掉相对丰度低于0.01%的物种。这并不代表Kraken2有严重缺陷任何基于k-mer的算法都难以区分共享高度保守序列的近缘物种你要做的是在灵敏度和特异度之间找到适合你项目的平衡点。6.5 数据库更新后结果无法复现Kraken2会定期更新预构建数据库不同时间点下载的库版本不一样里面参考基因组的版本、数量和taxonomy信息都会有变化。同一样本在不同版本数据库下跑出的结果不同这属于正常现象但也意味着做正式分析时必须冻结数据库版本。我一般会在项目目录下放一个environment.yaml文件里面同时记录Kraken2版本和数据库的下载日期、URL或MD5值。有的朋友会在论文method部分直接写Kraken2 v2.1.3 against Standard database downloaded on 2024-03-15这是业内比较认可的做法。7. 把Kraken2Bracken接入完整分析流程几个进阶思路7.1 从分类结果中提取目标reads做后续分析Kraken2运行时的--classified-out参数保留了所有分类成功的reads这让按物种筛选reads变得异常简单。比如我想提取样本中所有的真菌reads做单独组装只需要grep -E ^ classified_1.fq | head这种基于文件名的粗筛不符合需求更好的做法是写一个脚本从kraken2.output里提取指定taxid对应的read名字再从classified_out文件里用seqkit过滤。只要跑一次Kraken2就能拿到一个按物种索引的reads库之后不管是组装特定物种的基因组还是统计特定基因的存在与否都特别方便。这也是我推荐大家优先处理原始reads而不是处理contig的原因。7.2 批量处理多个样本时的流程自动化如果你有几十个样本一条条手动跑Kraken2不现实写个shell循环是最简单的for R1 in *_R1.fastq.gz; do R2${R1/_R1.fastq.gz/_R2.fastq.gz} SAMPLE${R1%%_R1.fastq.gz} kraken2 --db krakendb \ --paired \ --threads 16 \ --confidence 0.2 \ --output ${SAMPLE}.kraken.out \ --report ${SAMPLE}.kraken.report \ $R1 $R2 bracken -d krakendb \ -i ${SAMPLE}.kraken.report \ -o ${SAMPLE}.bracken \ -r 150 -l S -t 16 done样本量再大一些比如上百个我就建议用Snakemake了因为失败重跑、断点续跑、并行资源控制这些事情shell脚本维护成本太高。Snakemake的rule定义也很直观Kraken2和Bracken各自的输入输出界定清楚加上--restart-times参数跑自动化流程会从容很多。7.3 结果可视化与下游统计Bracken的输出文件接各种工具都比较自然。我做相对丰度柱状图时习惯把结果整理成长表用R的ggplot2画堆叠柱状图。做差异分析时可以把Bracken的fraction_total_reads矩阵作为输入接LEfSe或者是DESeq2性状分析。Kraken2的report文件也可以直接喂给Pavian这类网页工具做交互式探索Pavian的桑基图对于从界到种的分类层级展示挺直观适合给合作者快速看结果。7.4 要不要同时试一下Centrifuge和KaijuKraken2不是唯一的快速分类方案Centrifuge和Kaiju也有各自的拥趸。我自己的横向对比经验是Centrifuge的内存占用低、分类速度快但它的输出格式和Kraken2差异比较大Kaiju对蛋白水平的分类更敏感适合远源同源检测。但从生态完善度来看Kraken2的数据库体系、下游工具链、Bracken的配套、还有社区教程数量明显比另外两个更成熟。对刚开始接触宏基因组物种注释的人我还是建议先吃透Kraken2Bracken这一条流程后续如果有特殊需求再横向扩展。数据库的预训练和定期更新是个体力活但一旦把这套流程沉淀下来后面的分析效率真的能高出一大截。我给课题组搭完这套流程之后同事们从一个样本等一晚上变成了中午提交下午就能看到物种组成。这种效率上的变化对一个生信支撑人员来说比什么花哨的算法都实在。
返回列表