Skip to content

基因组注释

直接将数据比对到参考的基因组库中，并注释出来。为了后续的可视化，使用kraken2的结果模式。运行过程非常吃内存。

centrifuge

使用centrifuge来进行物种注释。

1
2
3
4
5
6
centrifuge \
    -x hpvc \
    -1 removeHuman/SRR10903401_1.fastq.gz \
    -2 removeHuman/SRR10903401_2.fastq.gz \
    -S centrifuge/SRR10903401.txt -p 8 \
    --report-file centrifuge/SRR10903401.report.tsv

SRR10903401.txt结果一共分为8列。从左到右，原始read ID，比对到数据库的序列ID（如果使用RefSeq或nt库则时AccessionID），物种分类ID，classification得分（序列之和），第二比对结果得分，比对序列长度，比对reads长度，reads比对上的物种数。

我们可以关注SRR10903401.report.tsv结果，结果中对物种进行了注释。

将结果转为kraken2的结果回报格式。

1
2
centrifuge-kreport -x hpvc centrifuge/SRR10903401.report.tsv \
    > centrifuge/SRR10903401.kreport.tsv

运行时报import imp Warning，因为python3.4后不再使用imp。需要把import imp改为import importlib。

kraken2

使用kraken2来进行物种注释。

1
2
3
4
5
6
7
kraken2 --db minikraken_8GB_20200312 \
    --threads 8 \
    --report kraken2/SRR10903401.kraken2.tsv \
    --output kraken2/SRR10903401.report.tsv \
    --gzip-compressed --paired \
    removeHuman/SRR10903401_1.fastq.gz \
    removeHuman/SRR10903401_2.fastq.gz

由于kraken2本次分析使用的数据库大小远比centrifuge的小，因此速度要快很多。

kraken2分析的下游还可以使用bracken来进行校正。bracken和kraken2可使用相同的数据库。

1
2
bracken -d minikraken_8GB_20200312 -i kraken2/SRR10903401.kraken2.tsv \
    -o kraken2/SRR10903401.bracken.txt -l S

其他

其他可用软件：

metaphlan3、Prokka、SURPI

思考

在获得结果后，一般来说在病原微生物的项目中，我们需要使用人源等浓度核酸同步实验作为一个阴性对照，再在分析结果中，除去阴性对照的结果。由于阴性对照的结果与待测样本结果的reads数数量级可能不同，因此很多时候可能要考虑RPM作为对照分析的值。另外，在实验时使用PCR-free的方式可能会对下游分析更好（但PCR-free会导致样本量极少，一般下限是50pmol）。在获得过滤结果后，再结合患者的临床表现进行分析，得到最终结果。

部分结果的可能会有同源性（包括与人类基因组的同源性），所以要报哪一条？在回报reads数时，建议回报unique map的reads。

其实mNGS我个人感觉对实验技术的依赖比下游生信分析更高。