博文

Kraken2 Vs qiime2 16S物种注释

已有 3328 次阅读 2021-2-16 17:55 |个人分类:biology|系统分类:科研笔记

最早接触Kraken2这个软件是在宏基因组，但官网上说其实这个软件也是可以用于16S物种注释的。当时没怎么在意，后面发现有个美国肠道微生物检测公司Thryve是使用这个软件进行物种注释的。最近发现2020年9月的一篇文章是比较了kraken2和qiime2的物种注释结果，详细见宏基因组公众号的文章。
在这里插入图片描述
今天我要说的是亲自拿一个肠道微生物的样本测试下这个软件的分类效果，再和自己用的qiime2流程获得的结果比较下，绝知此事要躬行呀。
首先说一个令我震惊的发现，作者发现使用qiime2进行物种分类的时候，使用1-4核心与6-32核心的结果不一样，也没发现具体原因。好吧，虽然基本使用8核心进行物种注释，但是还是觉得如果这样都影响结果，软件的稳定性还是值得警惕的。那么后面，我们应该尽量使用>6核心进行分析，防止可重复性受到这个问题的影响。

Kraken2的过程

# 激活环境source  ~/miniconda3/bin/activate
conda activate kraken2# 下载数据库wget ftp://ftp.ccb.jhu.edu/pub/data/kraken2_dbs/16S_Greengenes13.5_20200326.tgz# 解压
 tar zxvf 16S_Greengenes13.5_20200326.tgz# 物种分类DBNAME=16S_Greengenes_k2db/
SAMPLE="s5630696"THREADS=8

kraken2 \
--db ${DBNAME} \
--threads ${THREADS} \
--report ${SAMPLE}.kreport2 \
--paired ${SAMPLE}_1.fastq ${SAMPLE}_2.fastq > ${SAMPLE}.kraken2# 获得相对丰度bracken -d ${DBNAME} -r 150 -l G -i ${SAMPLE}.kreport2 -o ${SAMPLE}_g.bracken1234567891011121314151617181920

qiime2的过程

# 这里只从物种注释开始，前面的见官方文档或者宏基因组公众号的中文翻译qiime feature-classifier classify-sklearn
--i-classifier classifier_gg13_8.qza
-i-reads rep_seqs_${SAMPLE}.qza
--o-classification ${SAMPLE}_classsified_gg.qza
--verbose
--p-n-jobs ${THREADS}# 导出门和属的结果for i in 2 6 
do qiime taxa collapse \
  --i-table table_final.qza \
  --i-taxonomy taxonomy.qza \
  --p-level $i \
  --o-collapsed-table table-l$i.qzadone#导出 tsvfor file in ./table-l*.qzado
  base=$(basename $file .qza)
  echo $base
  qiime tools export \
  --input-path $file \
  --output-path $base
  biom convert --to-tsv -i $base/feature-table.biom -o $base.tsvdone123456789101112131415161718192021222324252627

比较结果

在门级别上，和文章中类似，基本没有差别，毕竟这个分类级别很高。
在这里插入图片描述
属水平上，差别就比较明显了，这意味着可能要选定一个方法一直用下去，对于微生物群落研究来说，保存、提取、建库以及16S的数据分析都会有所影响，所以每个研究的结果只能对每个研究负责，所以矛盾的研究结果也是不足为奇吧！结果还是基本一致的。