[bioinfo]kraken2を使ってコンタミリードを排除する
超並列シークエンサーのデータにおいてコンタミに由来するリードは解析をする際に様々な問題を引き起こす。de novoアセンブルをする際にも数パーセントは目的以外のリードが含まれるため、これを取り除くことを試みた。
kraken2とkraken-toolsのインストール。kraken2はjellyfishも必要とするがこれは別にpathが通った場所にインストール済みとする(今回はcondaでインストールした)。
git clone https://github.com/DerrickWood/kraken2 cd kraken2 ./install_kraken2 cd ../ git clone https://github.com/jenniferlu717/KrakenTools cd KrakenTools chmod +x *.py
kraken2-buildでkraken2のデータベースを構築する。kraken2にはカスタムで自分で決めたデータやncbiに登録されているデータなどをデータベースに登録できるが、今回はヒトゲノムと微生物ゲノムのコンタミを検出したかったので、standardを使った。データ量がかなり大きいのでこれだけでも半日かかる。
kraken2-build --standard -db krakendb --threads 8
データベースができたらkraken2を実行する。
kraken2 --threads 12 --quick --db krakendb --paired read_1.fq read_2.fq --report kraken2_output_report > kraken2_out
kraken-toolsのextract_kraken_reads.pyで元のデータからコンタミを取り除いたものを抽出する。今回はhumanでも、細菌、ウィルスでもないリードを取り除きたかったので、classifiedされていない、unknownのリードを抽出した。–maxを指定しないと100 Mリード以上はオミットされてしまうようだ。
extract_kraken_reads.py -k kraken_out -s1 read_1.fq -s2 read_2.fq -o read_krakened_1.fq -o2 read_krakened_2.fq -t 0 --fastq-output --max 300000000
2023/03/30 追記
nanoporeを使った16Sシーケンスのためにsilvaの16Sデータベースを構築した。
kraken2-vuild -db silva_16s --type silva
しかしながら、standardでのclassificationの方が成績が良かった(例えば、あるサンプルはstandardで80% classifiedに対して、16S_silvaは24%程度であった)。