kraken2就不介紹了，是一款挺好用的快速比對宏基因組軟件。
現(xiàn)在有了nt數(shù)據(jù)庫下面animal的序列，那么如何構(gòu)建kraken2的數(shù)據(jù)庫呢？
首先推薦把kraken2安裝到單獨(dú)的conda環(huán)境中，而不是把kraken2直接安裝
conda install kraken2 ? ? ? ? ? ? ? ? #不推薦
這種雖然方便，但是很可能與其它軟件產(chǎn)生環(huán)境變量沖突
看了看kraken2的說明書，它提供”archaea”, “bacteria”, “plasmid”, “viral”, “human”, “fungi”, “plant”, “protozoa”, “nr”, “nt”, “env_nr”, “env_nt”, “UniVec”, “UniVec_Core”數(shù)據(jù)庫，但確實(shí)沒有動(dòng)物庫
那么自己構(gòu)建一個(gè)庫吧。
1 首先利用kraken2自身下載數(shù)據(jù)庫的功能，下載taxonomy分類庫
kraken2-build --download-taxonomy --threads 24 --db $DBNAME
2 添加nt.animal.fa
kraken2-build --add-to-library nt.animal.fa --db $DBNAME
這一步就會(huì)報(bào)錯(cuò)，因?yàn)閚t.animal.fa庫里有一些奇奇怪怪的accession，例如4W1Z_7
其實(shí)這些accession也沒問題，但kraken2就是識(shí)別不出來，好吧，那只能更改源代碼
找到kraken2安裝的目錄，修改一下libexec/scan_fasta_file.pl
注釋掉die那一行，就OK了，吐槽一下kraken2編程真嚴(yán)謹(jǐn)，如果是我這句估計(jì)都不寫
另外，需要注意的是，nt.animal.fa庫里有些很短的序列，可能只有幾十bp,也可以寫個(gè)perl把這些序列刪了，我試了試，152G的文件變成了149G，似乎用處不大
3 構(gòu)建個(gè)性化庫
kraken2-build --build --db?$DBNAME ? ? ? ? ? #當(dāng)時(shí)自己好像沒加線程，以后試試加線程
4 運(yùn)行kraken2
kraken2 --db $DB --threads 20 --paired ../s-592_NDSW56977_1.fq ../s-592_NDSW56977_2.fq --report report --output output --classified-out cseqs#.fq --unclassified-out de-s-592_NDSW56977#.fq ? ? ? ? ? ?#kraken2好處就是還能把未分類的序列paired，真優(yōu)秀，這些序列就是去除宿主的序列了
