macでインフォマティクス

macでインフォマティクス

HTS (NGS) 関連のインフォマティクス情報についてまとめています。

複数のbiningツールを統合し、包括的なメタゲノム解析を行うパイプライン metaWRAP

2018 タイトル修正, 説明追加, step5エラー修正

2019 データベース作成ケアレスミス修正, ヒートマップ追加, インストールコマンド修正

2020 gzip圧縮fastqは使えないことを追記, インストール手順修正, salmonの定量ステップ修正

2021 4/27,4/28 dockerリンク追記

2021 5/1 my docker imageを修正

2021 6/16 binningのコード追加

 

 全メタゲノム(WMG)ショットガンシーケンシングによる微生物群集の研究は、それらの分類学的組成に加えて、微生物の代謝ポテンシャルの研究のための新しい道を開くものである[論文より ref.1-3]。これは、ヒトの健康、廃棄物処理、農業、および環境管理に応用され、機能的相互作用、抗生物質耐性、および微生物の集団動態を解釈し予測する能力を大幅に改善する[ref.4-6]。数百から数千のコミュニティメンバーからなるWMGショットガンシーケンシングリードは、データ分析と解釈のためのユニークな課題を生成する[ref,3、7]。 WMGデータ分析用ソフトウェアは、その数や複雑さが増し、そのようなデータを処理、分析、解釈する能力が向上している[ref.8-12]。しかし、これらのツールは、生物学者が作業するには厄介になる。 WMGの分野が拡大するにつれて、メタゲノムデータの統一された分析が行える包括的かつアクセス可能なソフトウェアが必要とされている[ref.7、11]。

WMG解析を実行するには、現在利用可能なツールを見つけて設定し、インストールして構成し、競合するライブラリと環境変数を指定し、次のツールに入力するための正しい形式に変換するスクリプトを作成する必要がある[ref.13,14] 。これらの課題は、メタゲノム分析を試みている人、特に計算経験のない研究者にとって、微生物ゲノミクスの進展を妨げている分野にとって大きな負担となる[ref.15]。既存の自動化されたパイプラインとクラウドサービスはモジュール性がなく、ユーザーが分析を制御できず、ゲノム解読メタゲノミクス、メタゲノムアセンブリのビニングによる推定されたゲノムの抽出に欠けていることが多い(ref.14,16-19)。

(一部略)ビニングに取り組むために、CONCOCT、MaxBin、およびmetaBATなどのいくつかの洗練されたツールが開発されているが、それはまだ活発に改善されている[ref. 9,19-21]。メタゲノムビンの品質で考えるべきことは、(1)completion、すなわち集団ゲノムのカバレッジレベル、および(2)contamination、すなわちこの集団に属していない別のゲノム配列の量、である。これらの測定基準は、各ビン内の普遍的なシングルコピー遺伝子を数えることによって推定できる[ref.22,23]。 CheckMは、ビンのtaxonomyのゲノムが持つと期待されるシングルコピー遺伝子を調べることで改善を試みる[ref.24]。ビン中に見出が予想されるシングルコピー遺伝子のパーセンテージがcompletionとして解釈される、重複して見出されるシングルコピー遺伝子のパーセンテージからcontaminationが推定される。

 ほとんどのメタゲノムビニングツールは、類似の配列特性、例えばK-mer組成およびコドン使用率、および複数のサンプルにわたる同様のリードカバレッジを有するscaffoldsを一緒にクラスタリングすることによってビンを抽出する[ref.25,26]。どのような場合でも単一のビニング方法が優れているわけではないため、ビン統合ツールは強みを組み合わせて、さまざまなアプローチの弱点を最小限に抑えようとする。 DAS_Toolは、提供されたすべてのビンセット内のシングルコピー遺伝子を予測し、異なるビニング予測からのビンを集計し、各集約からより完全なコンセンサスビンを抽出して、結果として得られるビンが最もduplicationが少なく、もっともシングルコピー遺伝子が多くなるようにする[ref.27](紹介)。この折り畳み方法は、ビンの完成を大幅に改善する。一方、Binning_refinerは、コンティグをより多くのビンに分割し、オリジナルのビンセットと比較して、それらが異なるビンにあるものが見つかれば、2つのコンティグが同じビン内に存在しないようにする(紹介)。これはコンティグをより多くのビンに分割し、contaminationを減らす[ref.28]。これらのアプローチはいずれも、異なる方法のビンセットを統合し、優れたビンセットをもたらすが、DAS_Toolはcontaminationを犠牲にcontaminationを増加させ、Binning_refinerは純度を優先させcontaminationは無くす。未知のドラフトゲノム品質を相対的に改善する別の方法は、所定のビンに属するリードを抽出し、メタゲノムの残りの部分とは別にそれらをアセンブリすることである。適切なベンチマークによって、このアプローチは、微生物群集の少なくともいくつかのビンの品質および下流アノテーションを有意に改善することができる。

 ショットガンメタゲノミクスの分野は比較的新しいため、メタゲノムビンを検査、分析、視覚化するソフトウェアが不足している。metagenomic scaffoldsのtaxonomyを正確に予測できるツールが存在するが(Taxator-tkなど)、メタゲノムビン全体を分類するツールはない[ref.29,30]。同様に、リードアライメントのでデプスに基づいてscaffoldsのカバレッジを評価する方法は多数あるが、多くのサンプルのビン全体のカバレッジを見つける方法はない[ref.31,32]。最後に、メタゲノミクス共同体の文脈でドラフトゲノムを視覚化するツールはない。 WMGデータ分析のための使いやすい統合されたツールの必要性、およびメタゲノムビン分析のための利用可能なツールの欠如が、MetaWRAP構築の動機を著者らに与えた。

 

 

f:id:kazumaxneo:20180927203210p:plain

f:id:kazumaxneo:20180927203346p:plain

MetaWRAPのワークフロー(上)およびbinning性能の比較(下)。このツールだけで、rawリードのQT、de novoアセンブリ、binning、binning結果の比較と統合、binned fasta由来リードを使っての再アセンブリ(パフォーマンスが上がる)、taxonomic assignment、prokkaによるアノテーションができるようになっている。付随して、分析をサポートする図も出力できる。ワークフロー図にあるように、目的によってステップを飛ばしたり、順番を変更したりすることが可能。図はGithubより転載。

 

下にmetaWRAPのモジュール一覧を載せた。様々なモジュールを備える(Githubより)。

f:id:kazumaxneo:20180927204018p:plain

 

 

metaWRAPに関するツイート


インストール

macのminiconda2.4.0でテストした。

依存のバージョン

https://github.com/bxlab/metaWRAP/blob/master/installation/dependancies.md

 

本体 Github

たくさんのツールを使うため、オーサーらのマニュアルでは、コンフリクトの心配がない仮想環境へ導入することが推奨されています。以下はREADMEの手順に沿っています。

#"metawrap-env"という名前で仮想環境を作成、pythonは2.7とする。mambaを使う。
mamba create -n metawrap-env python=2.7.14 -y
conda activate metawrap-env

#チャネルを以下の順で指定して追加
conda config --add channels defaults
conda config --add channels conda-forge
conda config --add channels bioconda
conda config --add channels ursky

#metawrapを追加
mamba install -c ursky metawrap-mg -y

2021 4/27 一部依存が対応しなくなってきている。オーサーが記載しているdocker imageを使う。(blastだけ動作しない。NCBIのデータベースが新しくなっているが、こちらは古いblastがインストールされているため)
docker pull quay.io/biocontainers/metawrap:1.2--0

#2021 4/29 my docker image (*3)
docker pull kazumax/metawrap:1.22
docker run -itv $PWD:/data --rm kazumax/metawrap:1.22
> conda activate metawrap-env
> metaWRAP -h

metawrap -h

Usage: metaWRAP [module]

Options:

 

read_qc Raw read QC module (read trimming and contamination removal)

assembly Assembly module (metagenomic assembly)

kraken KRAKEN module (taxonomy annotation of reads and assemblies)

blobology Blobology module (GC vs Abund plots of contigs and bins)

 

binning Binning module (metabat, maxbin, or concoct)

bin_refinement Refinement of bins from binning module

reassemble_bins Reassemble bins using metagenomic reads

quant_bins Quantify the abundance of each bin across samples

classify_bins Assign taxonomy to genomic bins

annotate_bins Functional annotation of draft genomes

 

--help | -h show this help message

--version | -v show metaWRAP version

--show-config show where the metawrap configuration files are stored

 

metawrap read_qc

$ metawrap read_qc

metawrap read_qc

 

Usage: metaWRAP read_qc [options] -1 reads_1.fastq -2 reads_2.fastq -o output_dir

Options:

 

-1 STR          forward fastq reads

-2 STR          reverse fastq reads

-o STR          output directory

-t INT          number of threads (default=1)

 

--skip-bmtagger dont remove human sequences with bmtagger

--skip-trimming dont trim sequences with trimgalore

--skip-pre-qc-report dont make FastQC report of input sequences

--skip-post-qc-report dont make FastQC report of final sequences

 

 

> metawrap assembly

$ metawrap assembly

metawrap assembly

 

Usage: metaWRAP assembly [options] -1 reads_1.fastq -2 reads_2.fastq -o output_dir

Options:

 

-1 STR          forward fastq reads

-2 STR          reverse fastq reads

-o STR          output directory

-m INT          memory in GB (default=24)

-t INT          number of threads (defualt=1)

-l INT minimum length of assembled contigs (default=1000)

 

--megahit assemble with megahit (default)

--metaspades assemble with metaspades instead of megahit (better results, but slower and required a lot of RAM)

 

 

> metawrap kraken 

$ metawrap kraken 

metawrap kraken

 

Run on any number of fasta assembly files and/or or paired-end reads.

Usage: metaWRAP kraken [options] -o output_dir assembly.fasta reads_1.fastq reads_2.fastq ...

Options:

 

-o STR          output directory

-t INT          number of threads

-s INT read subsampling number (default=all)

 

Note: you may pass any number of sequence files with the following extensions:

*.fa *.fasta (assumed to be assembly files) or *_1.fastq and *_2.fastq (assumed to be paired)

 

> metawrap binning

$ metawrap binning

metawrap binning

 

------------------------------------------------------------------------------------------------------------------------

-----                             Non-optional parameters -a and/or -o were not entered                            -----

------------------------------------------------------------------------------------------------------------------------

 

 

Usage: metaWRAP binning [options] -a assembly.fa -o output_dir readsA_1.fastq readsA_2.fastq ... [readsX_1.fastq readsX_2.fastq]

Note: you must chose at least one binning method. Or all at once!

Options:

 

-a STR          metagenomic assembly file

-o STR          output directory

-t INT          number of threads (default=1)

-m INT amount of RAM available (default=4)

 

--metabat2      bin contigs with metaBAT2

--metabat1 bin contigs with the original metaBAT

--maxbin2 bin contigs with MaxBin2

--concoct bin contigs with CONCOCT (warning: this one is slow...)

--run-checkm immediately run CheckM on the bin results (required 40GB+ of memory)

--single-end non-paired reads mode (provide *.fastq files)

 

> metawrap bin_refinement 

$ metawrap bin_refinement 

metawrap bin_refinement

 

------------------------------------------------------------------------------------------------------------------------

-----                             Non-optional parameters -o and/or -A were not entered                            -----

------------------------------------------------------------------------------------------------------------------------

 

 

Usage: metaWRAP bin_refinement [options] -o output_dir -A bin_folderA [-B bin_folderB -C bin_folderC]

Note: the contig names in different bin folders must be consistant (must come from the same assembly).

 

Options:

 

-o STR          output directory

-t INT          number of threads (default=1)

-m INT memory available (default=40)

 

-A STR folder with metagenomic bins

-B STR another folder with metagenomic bins

-C STR another folder with metagenomic bins

-c INT minimum % completion of bins [should be >50%] (default=70)

-x INT maximum % contamination of bins that is acceptable (default=10)

 

--skip-refinement dont use binning_refiner to come up with refined bins based on combinations of binner outputs

--skip-checkm dont run CheckM to assess bins

--skip-consolidation choose the best version of each bin from all bin refinement iteration

--keep-ambiguous for contigs that end up in more than one bin, keep them in all bins (default: keeps them only in the best bin)

--remove-ambiguous for contigs that end up in more than one bin, remove them in all bins (default: keeps them only in the best bin)

 

データベースの準備
目的に合わせてデータベースをダウンロードする。例えばmetaWRAPを使ってヒトゲノムのコンタミを除くread_qcを行いたいなら、hg38のデータベースをダウンロードする必要がある。

f:id:kazumaxneo:20180927204903p:plain

Githubより

 

全体のフローを確認したいので、今回は全データベースをダウンロードする(link)。

#checkmデータベース
mkdir MY_CHECKM_FOLDER
checkm data setRoot #プロンプトが出てMY_CHECKM_FOLDERのフルパスを聞かれるので回答する
# Now manually download the database:
cd MY_CHECKM_FOLDER
wget https://data.ace.uq.edu.au/public/CheckM_databases/checkm_data_2015_01_16.tar.gz
tar -xvf *.tar.gz
rm *.gz

#krakenデータベース
mkdir MY_KRAKEN_DATABASE
kraken-build --standard --threads 24 --db MY_KRAKEN_DATABASE
kraken-build --db MY_KRAKEN_DATABASE --clean

#NCBI nt BLASTデータベース
mkdir NCBI_nt
cd NCBI_nt
wget "ftp://ftp.ncbi.nlm.nih.gov/blast/db/nt.*.tar.gz"
for a in nt.*.tar.gz; do tar xzf $a; done

#NCBI taxonomyデータベース
mkdir NCBI_tax
cd NCBI_tax
wget ftp://ftp.ncbi.nlm.nih.gov/pub/taxonomy/taxdump.tar.gz
tar -xvf taxdump.tar.gz

#humanデータベース
mkdir BMTAGGER_INDEX
cd BMTAGGER_INDEX
wget ftp://hgdownload.soe.ucsc.edu/goldenPath/hg38/chromosomes/*fa.gz
gunzip *fa.gz
cat *fa > hg38.fa
rm chr*.fa
#index(かなりメモリーが必要。著者らは-Mで100GB与えている)
bmtool -d hg38.fa -o hg38.bitmask
srprism mkindex -i hg38.fa -o hg38.srprism -M 100000

 

データベースダウンロード後、config-metawrapのパスを修正する。

#configのパスを調べる
which config-metawrap

 
#場所がわかったらemacsなどで開いてデータベースパスを修正
emacs /home/kazu/.pyenv/versions/miniconda3-4.0.5/envs/metawrap-env/bin/config-metawrap

#例えばKrakenなら
KRAKEN_DB=/path/to/my/database/MY_KRAKEN_DATABASE 

#同様に
BLASTDB=/your/location/of/database/NCBI_nt
TAXDUMP=/your/location/of/database/NCBI_tax

詳細リンク

 

テストラン

https://github.com/bxlab/metaWRAP/blob/master/Usage_tutorial.md

 

準備

EMBL-EBIのシーケンスリードアーカイブEuropean Nucleotide Archiveからシーケンシングデータをダウンロードし、解凍する。

wget ftp.sra.ebi.ac.uk/vol1/fastq/ERR011/ERR011347/ERR011347_1.fastq.gz 
wget ftp.sra.ebi.ac.uk/vol1/fastq/ERR011/ERR011347/ERR011347_2.fastq.gz

wget ftp.sra.ebi.ac.uk/vol1/fastq/ERR011/ERR011348/ERR011348_1.fastq.gz
wget ftp.sra.ebi.ac.uk/vol1/fastq/ERR011/ERR011348/ERR011348_2.fastq.gz

wget ftp.sra.ebi.ac.uk/vol1/fastq/ERR011/ERR011349/ERR011349_1.fastq.gz
wget ftp.sra.ebi.ac.uk/vol1/fastq/ERR011/ERR011349/ERR011349_2.fastq.gz

gunzip *qz
mkdir RAW_READS
mv *fastq RAW_READS

こうなっているはず。

f:id:kazumaxneo:20181002102129j:plain

 

1) Read_QC: read trimming and human read removal

クオリティトリミング(QT)を実行し、QT前後のクオリティレポートを出力する。

mkdir READ_QC

metawrap read_qc -1 RAW_READS/ERR011347_1.fastq -2 RAW_READS/ERR011347_2.fastq
-t 40 -o READ_QC/ERR011347

metawrap read_qc -1 RAW_READS/ERR011348_1.fastq -2 RAW_READS/ERR011348_2.fastq
-t 40 -o READ_QC/ERR011348

metawrap read_qc -1 RAW_READS/ERR011349_1.fastq -2 RAW_READS/ERR011349_2.fastq
-t 40 -o READ_QC/ERR011349
  • -1    forward fastq reads

  • -2    reverse fastq reads

  • -o    output directory

  • -t     number of threads (default=1)

この後の作業に備え、QT後のfastqを/CLEAN_READSに移動。

mkdir CLEAN_READS
for i in READ_QC/*; do
b=${i#*/}
mv ${i}/final_pure_reads_1.fastq CLEAN_READS/${b}_1.fastq
mv ${i}/final_pure_reads_2.fastq CLEAN_READS/${b}_2.fastq
done

注意: gzip圧縮されたfastqは使えない。またペアエンドはxxx_1.fastq、xxx_2.fastqなどになっている必要がある。

 

2) Assembly: metagenomic assembly and QC with metaSPAdes or MegaHit

巨大なデータでない限り、De novoアセンブリにはmegahitよりmetaspadesが推奨されている。

cat CLEAN_READS/ERR*_1.fastq > CLEAN_READS/ALL_READS_1.fastq
cat CLEAN_READS/ERR*_2.fastq > CLEAN_READS/ALL_READS_2.fastq

#スレッドは24、メモリは100GB指定
metawrap assembly -1 CLEAN_READS/ALL_READS_1.fastq -2 \
CLEAN_READS/ALL_READS_2.fastq -m 100 -t 40 --metaspades -o ASSEMBLY
  •  -m   memory in GB (default=24)

 

3) Kraken: taxonomy profiling and visualization or reads or contigs

krakenを動かしtaxonomy profilingを実行する。

metawrap kraken -o KRAKEN -t 40 -s 1000000 CLEAN_READS/ERR*fastq ASSEMBLY/final_assembly.fasta
  • -s   read subsampling number (default=all) 

f:id:kazumaxneo:20181002223943j:plain

krakenのtaxonomic profiling 結果は、内部でkronaを使って可視化される。

f:id:kazumaxneo:20181002223837j:plain

 

 

4) Binning: initial bin extraction with MaxBin2, metaBAT2, and/or CONCOCT

3ツールのbinningを実行する(binning1回目)。(*2)マージしたデータと、マージせず個別に指定したデータで結果は変わる。一般に後者の方がbinned.fastaは増える。

metawrap binning -o INITIAL_BINNING --metabat2 --maxbin2 --concoct\
-t 40 -a ASSEMBLY/final_assembly.fasta\
CLEAN_READS/ERR*fastq

#concoctの代わりにmetabat1を使う
metawrap binning -o INITIAL_BINNING --metabat2 --maxbin2 --metabat1\
-t 40 -a ASSEMBLY/final_assembly.fasta\
CLEAN_READS/ERR*fastq

#ペアエンドリードを個別に指定、例えばsampleA、B、C、Dがあるとする
metawrap binning -o INITIAL_BINNING \
-t 40 --metabat2 --maxbin2 --concoct -a scaffolds.fasta \
sampleA_1.fastq sampleA_2.fastq sampleB_1.fastq sampleB_2.fastq \
sampleC_1.fastq sampleC_2.fastq sampleD_1.fastq sampleD_2.fastq

Unable to find proper fastq read pair in the format *_1.fastq and *_2.fastq 

というエラーが出たら、fastqの名前を末尾が *_1.fastq、 *_2.fastqになるよう修正して下さい。

  • metaspadesハイブリッドアセンブルして得たより長いsaffolds.fastaも使える(metaspadesをmetaWRAPを使わずにランする)。

 

このステップは特にリソース集約的な作業になる。マニュアルでは、データが巨大なら個別にランした方が良いかもしれないと述べられている。

出力ディレクト

f:id:kazumaxneo:20181002215438j:plain

 

 

5) Consolidate bin sets with the Bin_refinement module

3つのbinを合体させる。ここでは3結果の統合だが、4以上ある時は複数回走らせる。checkmを使う。まだcheckmの設定を行なっていないなら、先にデータベースの準備のcheckmのコマンドを実行しておく。

metawrap bin_refinement -o BIN_REFINEMENT -t 40 \
-A INITIAL_BINNING/metabat2_bins/ \
-B INITIAL_BINNING/maxbin2_bins/ \
-C INITIAL_BINNING/concoct_bins/ \
-c 50 -x 10
  • -c   minimum % completion of bins [should be >50%] (default=70)

  • -x   maximum % contamination of bins that is acceptable (default=10)

パラメータのminimum compleiton (-c)とmaximum contamination (-x) の最適値は目的によって変わり得る(The high completion but high contamination, or the less complete but more pure bin? (Githubより))。テストでは、カバレッジが低いため、-c 50 -x 10にしている。仮に-c 90 -x 5にすれば、コンタミは減りより完全に近いbinの割合が増えるが、bin数は減る。

出力ディレクト

f:id:kazumaxneo:20181013130647j:plain

✳︎注  指定のcheckm version1.0.7を使ったが、初期はエラーが出て修正できなかったため--skip-checkmをつけてcheckmを飛ばして走らせていた。その後、pplacerのメモリ不足が原因でcheckmが正常に動いていないことがわかった。dockerコンテナ立ち上げ時に80gほどメモリ指定し、CPUコア数も減らしてメモリ使用量を抑えることでエラーなくランできるようになった。

 

  

6) Visualize the community and the extracted bins with the Blobology module

contigを可視化する。

metawrap blobology -o BLOBOLOGY --bins BIN_REFINEMENT/metawrap_bins\
-t 40 -a ASSEMBLY/final_assembly.fasta \
CLEAN_READS/ERR*fastq

出力ディレクト

f:id:kazumaxneo:20181015182540j:plain

final_assembly.binned.blobplotの中身(excelで表示)

f:id:kazumaxneo:20181003144947j:plain

final_assembly.binned.blobplot.taxlevel_order.png

f:id:kazumaxneo:20181015180633p:plain

final_assembly.blobplot.taxlevel_phylum.png

f:id:kazumaxneo:20181015180756p:plain

final_assembly.blobplot.taxlevel_superkingdom.png

f:id:kazumaxneo:20181015180830p:plain

final_assembly.blobplot.bin.png

f:id:kazumaxneo:20181015180846p:plain

 

またはBlobsplorerを使って可視化する(github)。

git clone https://github.com/mojones/blobsplorer.git
cd blobsplorer/
open .

blobsplorer.htmlを開く。

f:id:kazumaxneo:20181003144130j:plain

Blobsplorerに得られたfinal_assembly.binned.blobplotを読み込ませる。

f:id:kazumaxneo:20181003144202j:plain

使用方法については右上のToggle Helpをクリック。

 

7) Find the abundaces of the draft genomes (bins) across the samples

salmonを使ってcontigを定量し、各binningゲノムの平均abundanceを計算する。5のBin_Refinment module結果を使うことが推奨されている。

metawrap quant_bins -b BIN_REFINEMENT/metawrap_bins -o QUANT_BINS \
-a ASSEMBLY/final_assembly.fasta sample1/*fastq sample2/*fastq sample3/*fastq

  "genome copies per million reads"で正規化されたabundance情報ファイルが出力される。

f:id:kazumaxneo:20181013130520j:plain

またabundance結果からクラスタリングして、量比をヒートマップで可視化したファイルも出力される。

f:id:kazumaxneo:20181013130545p:plain

追記

より大きなデータセット

f:id:kazumaxneo:20190415133739p:plain

 

8) Re-assemble the consolidated bin set with the Reassemble_bins module 

各bin由来リードを集め、個別にリアセンブリする。精度がより高まるとされる。bwaでリードをマッピングし(ミスマッチ0の厳密な条件と、ミスマッチ<5のエラー許容条件の2種類)、ペアエンドを回収(片方しかマッピングされないペアエンドも両方回収される)、spades(metaspadesより連続性改善が期待できる)で再アセンブリする。

metawrap reassemble_bins -o BIN_REASSEMBLY \
-t 40 -m 800 -c 50 -x 10 -b BIN_REFINEMENT/metawrap_bins \
-1 CLEAN_READS/ALL_READS_1.fastq -2 CLEAN_READS/ALL_READS_2.fastq

アセンブリ結果は出力されるBIN_REASSEMBLY/reassembled_bins.statsを確認する。厳密な条件のマッピングで得られたコンティグは名前にstrict、 ミスマッチ<5のエラー許容条件で得られたコンティグは名前にpermissive が付く。

 

出力ディレクト

f:id:kazumaxneo:20181013190128j:plain

またcheckMの結果をまとめたPNGファイルが出力される(IN_REASSEMBLY/reassembled_bins.png)。

reassembled_bins.png

f:id:kazumaxneo:20181013190159p:plain

reassembly_results.png

f:id:kazumaxneo:20181013190218p:plain

 

 

9) Determine the taxonomy of each bin with the Classify_bins module 

Taxator-tk(リンク)を使い、各contigにphylogenic情報をアサインする。ダウンロードしたNCBI_nt と NCBI_tax databases が使われている(精度は使用されるデータベースに依存する)。

metawrap classify_bins -b BIN_REASSEMBLY/reassembled_bins \
-o BIN_CLASSIFICATION -t 40

 分類結果はBIN_CLASSIFICATION/bin_taxonomy.tabに出力される。

出力ディレクト

f:id:kazumaxneo:20181015135536j:plain

 

10) Functionally annotate bins with the Annotate_bins module

アノテーションを行う。Prokkaが使われている。

metaWRAP annotate_bins -o FUNCT_ANNOT -t 40 \
-b BIN_REASSEMBLY/reassembled_bins/

Prokkaの最終出力はFUNCT_ANNOT/prokka_out に出力される。

出力ディレクト

f:id:kazumaxneo:20181015135751j:plain

FUNCT_ANNOT/prokka_out/ binの出力

f:id:kazumaxneo:20181015135903j:plain

 

 

しかしすごいツールですね。オーサーらのツイートのやり取りを見ると、biocondaへの登録を進めるようなので、そのうちbiocondaとdockerコンテナを介した利用もできるようになりそうです。

引用

MetaWRAP—a flexible pipeline for genome-resolved metagenomic data analysis

Gherman V. Uritskiy, Jocelyne DiRuggiero, James Taylor

Microbiome 2018 6:158

 

関連ツール


 

 

*1

checkmのエラーメッセージから 、メモリ要求量の問題か、一部のビニングツールのランが失敗しているか、checkmのデータベースが正しく構築できていないことがあたりに問題があるかと考えた。対処療法としていくつかサンプルを変えて試してみたが解決できなかった。最終的にはメモリ不足でpplacerがcore dumpしていることがわかった。→ メモリを増やすことでランできるようになった。

 

 *2

concoctのopen MPIエラーでループしてしまったら以下のスレッドのアドバイス通り、

https://github.com/BinPro/CONCOCT/issues/232

metawrapの仮想環境で

#concoctをアンインストール

conda uninstall concoct

#concoctのコードを取ってきて

git clone https://github.com/BinPro/CONCOCT.git

cd CONCOCT

#CONCOCTのインストール

conda install -c conda-forge openblas=0.3.3

python setup.py install

そのあとmetawrapを導入し直すことで直った。

conda install -c ursky metawrap-mg -y

 

 

*3 

blastコマンドとNCBI のBLASTデータベースのバージョンが異なると、つまりインストールしているblastコマンドと異なるバージョンのblastコマンドで作成されたデータベースをダウンロードしてしまうと、データベースを認識できずblastのステップでエラーになる。

最新のBLASTデータベースをダウンロードしたなら対策は簡単で、別の仮想環境にblastコマンドの最新版をダウンロードして、現在のblastnパスに上書きする形でそのバイナリのパスを張れば解決する。すなわち、

blastを一度アンインストールする。metawrapのパスも除外されるが、依存関係のパスは通ったままなので、再度metawrapを導入する。これでblastのみ最新版になる。もちろんソースからビルドして置き換えてもO.K。

conda activate metawrap-env

#blastのバージョンチェック
blastn -version

#blastを除く
mamba remove blast

#mambaでblast最新版を導入(bioconda)
mamba install -c bioconda blast -y
blastn -version

#metawrapを改めて導入(ダウンロード済みプログラムにパスが通るだけなので2度目は早い)
mamba install -c ursky metawrap-mg -y

#blastのバージョンチェック
blastn -version

上の流れで対策できる。ただし上で紹介しているイメージはこの方法で最新のblastを導入できなかったため、ソースコード(ncbi-blast-2.11.0+-src.tar.gz )からビルドしてcondaのパスのblastを置き換えている。

metawrap classify_binsのmegablastでエラーになる時は

conda remove blast
mamba install -c ursky metawrap-mg -y

 で、失敗した出力を指定して再びランすることで(途中から再開される)、一応最後までランすることができた。