macでインフォマティクス

macでインフォマティクス

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

メタゲノム解析ツール

 使ってみて便利だったツールを紹介する。

 

Genome sequences of rare, uncultured bacteria obtained by differential coverage binning of multiple metagenomes Albertsen et al. (2013)

 メタゲノムデータから、各生物ごとのデータを大まかに仕分け、その後の処理、例えばドラフトゲノム構築などを助ける手法。難培養性でデータ比率がシーケンス比率が全体の1%程度の株のゲノムでも、ドラフトゲノムを構築できるとされる。初めて聞いた時はシンプルかつエレガントなやり方に目から鱗が落ちた。

 

 ポイントは、ゲノムDNAを抽出比率が変化するような2つのやり方で抽出することにある。論文中では下水処理の活性化汚泥からフェノール有り/無しでDNAを抽出して、異なるindexをつけilluminaの150 bp x 2シーケンスしている。得られた1.7億リードと9千万リードをde novoアセンブルする。

 そこから先が面白い。全コンティグのカバレッジをフェノール有り/無しそれぞれで算出し横軸と縦軸にとってプロットする。さらに各コンティグのGC含量も計算して、各contigのプロットに色をつける。こうするとどうなるかというと、同じような色のコンティグが特定の領域に収束する(下の図)。

f:id:kazumaxneo:20170523104958j:plain

公式ページから引用。

グラフは4つの要素が表現されている。縦軸と横軸はフェノール有り/無しのカバレッジ、プロットの大きさはコンティグの長さ(bp)、そしてプロットの色はGG含量となる。Rのggplot2で作成されている。

 

 特定の領域に集まった同じような色のプロットは、同一生物か、極めて似た生物のゲノムのコンティグと考えられる。理由は、DNAの抽出率は細胞壁構造の違いなどにより生物ごとに一定の値になるためである。  

 例えば、あるDNA抽出条件では生物Aの膜は100%破れDNAが抽出されるが、生物Bの膜は50%しか破れないため、Aの半分のコピー数しかDNAが抽出されない(AとBが同一細胞数いると仮定して)。もう一つのDNA抽出条件では割合が逆転し、Aは50%破れ、Bは100%破れる。それらのシーケンスデータをde novoアセンブルして、できたコンティグそれぞれのカバレッジを計算して縦軸横軸にプロットするわけである。

 グラフができたら、グラフの中で1箇所に集まった同じ色のコンティグをラフに囲み回収する(ビニングする)。回収後、コンティグをさらに伸ばすために、再アセンブルしたりコンティグ間にまたがってマップされるペアリードを探す。そこから元にコンティグの隣接関係を絞り込み、edge-nodeのグラフとして可視化する。

 可視化することでedgeが切れている部位があれば、別生物のデータが混じっているか、情報が足りないか(binning時にロスしたか、シーケンス不足)素早く確認することができる。その菌のゲノムを決めたければ、node間をPCRしたり、追加でロングリードシーケンスしてフィニッシングまで持っていくことも可能となる。  

 

 手順はこちらに非常に丁寧に書かれている。手法についてのスクリプトやRのコマンドについて丁寧に教えてくれるのでインフォマティクスの勉強にもなる。

 

 

 

PhyloPhlAn Segata et al. (2014)

 

  比較する生物間に保存されている全orfを使って系統樹を作成する手法。16Sのような単独のorfを使った系統樹では、比較対象が近縁種すぎると距離がゼロになる可能性がある。この手法は保存されたorfを全て使い、より精度の高い分子系統樹を自動で書いてしまうというものである。

 保存されている全遺伝子を使うため、わずかな遺伝子の違いでも評価できる。株間の系統解析のような,同一種の分析にも有用と思われる。 

 

解析の流れはこちらに書かれている。中のリンクからダウンロードする。

ランにはFasttreeが必要。公式ページこちらからダウンロードする。またはbrewでインストールしてもよい。

brew install FastTree

 

アミノ酸ファイル(faa)をダウンロードして一箇所のフォルダに集め、そのフォルダ名をinputとして実行する。

./phylophlan.py -u example --nproc 12

  -u, --user_tree: Build a phylogenetic tree using user genomes only 

  --nproc N:  The number of CPUs to use for parallelizing the blasting

 

 

exampleデータを実行するとoutput/exampleにnwk形式のファイルが出力される。系統樹ソフトMEGAに読み込ませて結果を確認。

f:id:kazumaxneo:20170523104808j:plain

 

 

 

ConStrains Luo et al. (2015)

メタゲノムデータから、遺伝子間のSNPなどを指標に種内間の多様性を評価するためのツール。ある種の生物がどのくらいの多様性を持ってどのような比率で存在しているか解析することが可能。

 

gitからダウンロードする。 (Wikiリンク)

git clone https://bitbucket.org/luo-chengwei/constrains

ランにはpytohn2.7環境が必要。さらにpythonのBiopython、Numpy、Scipy、NetworkXを必要とする。そのほか、Bowtie、Metaphlan2などが必要。ないものだけ以下のコマンドでインストール。

pip install Biopython
pip install Numpy
pip install Scipy
pip install networkx
brew install Bowtie

Metaphlan2は2015年に発表された手法。ダウンロードリンクは下のMetaphlan2紹介の方に記載している。

 

テストコマンド

はじめに以下のコマンドでfastqデータとconfig_fileを自動ダウンロードする。

python utils/download_tutorial.py -o test #30分くらいかかる。

ダウンロードされたfastqはE.coliのいくつかの株の混ぜ物からARTでシミュレートされたペアリードデータらしい。データの比率は以下のように書いてある。

The two samples contain ART simulated E.coli Illumina reads. In sample_1, the composition is 5% (ETEC H10407) and 95% (HS), and in sample_2, the composition is 15.2% (IAI1), 26.3% (CFT073), and 58.5% (O157:H7).

テストではこの株の同定を試みている。

解凍する。 

tar xzvf constrains.test/fq.tar.gz

できた4つのfastqファイルをターミナルのカレントディレクトリに移動する。またはダウンロードされたtest.confを編集して、fqのパスを修正する。中身を見ると

>cat test/test.conf
// sample: sample_1
fq1: ./sample_1.1.fq
fq2: ./sample_1.2.fq
//
sample: sample_2
fq1: ./sample_2.1.fq
fq2: ./sample_2.2.fq

このようにfqファイルのパスとサンプル名が書かれている。

 

テストランは以下のように行う。

python ConStrains.py -c test/test.conf -o Constrains_test -t 12

 -c: コンフィグファイル

-t: スレッド数

 

テストは10分くらいで終わる。Constrains_test/results/Intra_sp_rel_ab.profilesを開くと

# Species   strain_ID   masked_samples  sample_1    sample_2

Escherichia_coli    str-1   NA  0.00000000  14.88058489

Escherichia_coli    str-2   NA  3.75299123  0.00000000

Escherichia_coli    str-3   NA  96.24700877 85.11941511

となっていた。Intra_sp_rel_ab.profilesは同一種内の株の比率 (%) で(つまり合計は常に100)、Overall_rel_ab.profilesはメタゲノムのコミュニティ全体から見た株の比率を表す。

 

SRAに登録されていったbiogass producing micorobialのメタゲノムデータでも試してみた。サンプル数1は12GBのfastqが2つある。-t=16でランすると1時間程度で終了した。結果のOverall_rel_ab.profilesを見ると下のようになっていた。

# Species strain_ID masked_samples ERR843249
Lactobacillus_curvatus NA,insufficient 1 0.28228
Enterococcus_cecorum NA,insufficient 1 2.77777
Coprococcus_catus NA,insufficient 1 0.11833
Weissella_koreensis NA,insufficient 1 0.37525
Prevotella_copri NA,insufficient 1 0.22341
Methanoculleus_bourgensis NA,insufficient 1 25.03193
Aminobacterium_colombiense NA,insufficient 1 1.07588
Streptococcus_macedonicus NA,insufficient 1 0.02517
Lactobacillus_johnsonii NA,insufficient 1 0.11143
Thermoanaerobacter_thermohydrosulfuricus NA,insufficient 1 5.57953
Streptococcus_lutetiensis NA,insufficient 1 0.6998
Lactobacillus_reuteri NA,insufficient 1 0.45737
Streptococcus_infantarius NA,insufficient 1 8.61296
Leuconostoc_citreum NA,insufficient 1 0.15594
Lactobacillus_amylovorus NA,insufficient 1 3.39372
Agrobacterium_sp_H13_3 NA,insufficient 1 34.99582
Aerococcus_viridans NA,insufficient 1 0.49274
Methanoculleus_marisnigri NA,insufficient 1 0.24412
Staphylococcus_lentus NA,insufficient 1 0.12412

 シーケンス率の異なるかなり多様な菌が混じっていることがわかる。

 

  

 

 

MetaPhlAn2 Luo et al. (2015)

 

 解説ページ。

Metaphlanはメタゲノムデータから種のレベルで分類する手法。MG-RASTやFOCUSと同様のツール。clade特異的なマーカー遺伝子を元に仕分けを行なっており、カバレッジを計算することでデータ中の比率も自動計算されるようになっている。Metaphlan2ではこれをさらに高精度化して、strainの分類まで検出することが可能になっている。7500種以上の生物から100万近いマーカー遺伝子を準備し、ウィルスも検出できるようになったらしい。 

こちらのリンクからダウンロードする(大量のマーカ遺伝子を含むため1GB以上のサイズになっている)。

ランにはpython2.7のライブラリが必要。

pip install numpy
pip install Biom-format

リードをアライメントさせるために、bowtie2にパスが通っている必要がある。

 

 

-hで表示されるコマンドを載せておく。

 

usage: metaphlan2.py --input_type

                     {fastq,fasta,multifasta,multifastq,bowtie2out,sam}

                     [--mpa_pkl MPA_PKL] [--bowtie2db METAPHLAN_BOWTIE2_DB]

                     [--bt2_ps BowTie2 presets] [--bowtie2_exe BOWTIE2_EXE]

                     [--bowtie2out FILE_NAME] [--no_map] [--tmp_dir]

                     [--tax_lev TAXONOMIC_LEVEL] [--min_cu_len]

                     [--min_alignment_len] [--ignore_viruses]

                     [--ignore_eukaryotes] [--ignore_bacteria]

                     [--ignore_archaea] [--stat_q]

                     [--ignore_markers IGNORE_MARKERS] [--avoid_disqm]

                     [--stat] [-t ANALYSIS TYPE] [--nreads NUMBER_OF_READS]

                     [--pres_th PRESENCE_THRESHOLD] [--clade] [--min_ab] [-h]

                     [-o output file] [--sample_id_key name]

                     [--sample_id value] [-s sam_output_file]

                     [--biom biom_output] [--mdelim mdelim] [--nproc N] [-v]

                     [INPUT_FILE] [OUTPUT_FILE]

metaphlan2.py: error: argument --input_type is required

user-no-MacBook-Pro:biobakery-metaphlan2-f1dcf3958459 user$

fastqを分析するなら

metaphlan2.py metagenome.fastq --bowtie2out metagenome.bowtie2.bz2 --nproc 12 --input_type fastq > profiled_metagenome.txt

 メタゲノムじゃないデータを使うと、以下のような単一の生物cladeを示したtextデータが出力される。正しく生物名を検出できている。

f:id:kazumaxneo:20170628151652j:plain

出力される数値はノーマライズされた100%相対値。存在率計算にこれ以上の計算は必要ない。出力についてはwiki参照。

  • --min_ab The minimum abundance (1.0).
  • -t <ANALYSIS TYPE>      Type of analysis to perform:
  1.  -t rel_ab profiling a metagenomes in terms of relative abundances
  2. -t rel_ab_w_read_stats profiling a metagenomes in terms of relative abundances and estimate the number of reads comming from each clade.
  3. -t reads_map mapping from reads to clades (only reads hitting a marker)
  4. -t clade_profiles normalized marker counts for clades with at least a non-null marker. To obtain all markers present for a specific clade and all its subclades, the '-t clade_specific_strain_tracker' should be used.
  5. -t marker_ab_table normalized marker counts (only when > 0.0 and normalized by metagenome size if --nreads is specified)
  6. -t marker_pres_table list of markers present in the sample (threshold at 1.0 if not differently specified with --pres_th [default 'rel_ab']

オプション一覧

ノーマライズされた値でなくマーカー遺伝子にヒットしたリード数を出力するなら、-t marker_countsを使う。

 

データの中にどんな生物のデータがどれだけあるのか見極めるには便利なツールではあるが、マーカー遺伝子とヒットしなかったリードは当然評価の対象外となる。例えば、非モデル生物でドラフトゲノムのデータも登録されていないような生物のリードが混じっているなら、許されないレベルの誤差が生じる可能性がある。論文化に当たってはqPCRの裏付けが必要かと思われる。

 

 

 

Metacoder Foster et al. (2017)

 Rのパッケージとして提供されている。

インストールは、Rを起動して

install.packages("metacoder")

  

制作途中

 

 

Mycc Lin et al. (2016)

ダウンロードリンク

 

 

 

 

 

 

 

------------------------------------------引用------------------------------------------

Genome sequences of rare, uncultured bacteria obtained by differential coverage binning of multiple metagenomes

Genome sequences of rare, uncultured bacteria obtained by differential coverage binning of multiple metagenomes. - PubMed - NCBI

 

ConStrains identifies microbial strains in metagenomic datasets.

ConStrains identifies microbial strains in metagenomic datasets. - PubMed - NCBI

 

PhyloPhlAn is a new method for improved phylogenetic and taxonomic placement of microbes

Nicola Segata, Daniela Börnigen, Xochitl C. Morgan, and Curtis Huttenhower. Nature Communications 4, 2013

https://www.ncbi.nlm.nih.gov/pubmed/23942190

  

MetaPhlAn2 for enhanced metagenomic taxonomic profiling

Duy Tin Truong, Eric A Franzosa, Timothy L Tickle, Matthias Scholz, George Weingart, Edoardo Pasolli, Adrian Tett, Curtis Huttenhower & Nicola Segata. Nature Methods 12, 902–903 (2015)

MetaPhlAn2 for enhanced metagenomic taxonomic profiling. - PubMed - NCBI

 

Metacoder: An R package for visualization and manipulation of community taxonomic diversity data

Metacoder: An R package for visualization and manipulation of community taxonomic diversity data

 

Accurate binning of metagenomic contigs via automated clustering sequences using information of genomic signatures and marker genes

Hsin-Hung Lin & Yu-Chieh Liao Scientific Reports 6, Article number: 24175 (2016) doi:10.1038/srep24175

https://www.nature.com/articles/srep24175