macでインフォマティクス

macでインフォマティクス

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

メタゲノムデータの平均ゲノムサイズや総カバレッジを推定する MicrobeCensus

metagenome evaluation tool library normalize

 

　ショットガンメタゲノミクスは、人体や環境の微生物群集の機能的構成を特徴づけるためにますます使用されてきている[論文より　ref.1-4]。これらの研究の共通の目標は、遺伝子ファミリー存在量を定量化し、環境、宿主の表現型、または実験条件の間で豊富さが異なる微生物の遺伝子またはパスウエイを同定することである。比較メタゲノミクスは、微生物が地球上の無数の環境にどのように機能的に適応しているか、そして微生物群集の機能的構成の変化がヒトの健康と疾病にどのように影響するかを明らかにしている[ref.5,6]。機能的変異は、同じショットガンデータから推定された分類学的変異の状況において解釈され、進化的および生態学的プロセスについての手がかりを与えている。これらの研究は、微生物群集DNAからの遺伝子量の正確な定量を必要とする。

　見過ごされがちだが、比較的メタゲノム解析を行う場合、微生物群における細胞の平均ゲノムサイズ（AGS）を考慮することが重要である。統計的な観点から、AGSは、コミュニティ間で豊富な遺伝子を比較する際に潜在的なバイアスの原因となる可能性があります。具体的には、コミュニティからある遺伝子をサンプリングする確率は、そのコミュニティのAGSに反比例する[ref.7 pubmed]。これを考慮しなければ、細胞あたり等コピー数で存在する遺伝子が異なるという結果を引き起こす（すなわち、偽陽性）。ライブラリーの大きさ（すなわちシーケンシングデプス）と遺伝子長に加えて、微生物群集間の生態学的に意味のあるゲノムの違いを同定するために、比較メタゲノム研究ではAGSを考慮する必要がある。 

　AGSは、環境内の微生物に作用する生態学的および進化的な力を理解するためにも重要である。生態学的な観点から、微生物のゲノムサイズは、環境の複雑さ、代謝生活様式、およびコミュニティのニッチを反映している可能性がある[ref.8-12]。例えば、腸内では、より大きいゲノムを有する生物はより一般的なライフスタイルに従うと考えられ、より小さなゲノムを有する生物はより特殊化されると考えられている[ref.8 pubmed]。例えばMetanobrevibacter smithii（1.9Mb）がH2の利用に特化したメタン菌であるのに対して、大きなゲノムを有するBacteroides thetaiotamicron（6.5Mb）はホストの食事およびホスト由来のグリカンの両方を利用できる代謝能を持っている。進化の観点から、ゲノムサイズの減少は、少数の集団における遺伝的浮動または大集団におけるgenome streamlining [ref.13,14]（引用　suggesting that reduction in genome size is beneficial and driven by positive selection）を反映する可能性がある。遺伝的浮動は、突然変異の欠失方向にバイアスが固定されることによって微生物ゲノムのサイズを形成することができる[15]。対照的に、genome streamliningは、栄養素が増殖をリミットするような大集団において、栄養素を効率的な使用するような選択圧がかかることでゲノムサイズの減少を促進する[ref.14]。

　現時点までにメタゲノミクスデータ用に設計されたソフトウェアがないため、微生物群のAGSを迅速かつ正確に推定することは不可能であった。これは、ヒト微生物を含む多くの環境におけるAGS変異の程度および影響の両方に対する我々の理解を制限している。唯一公的に利用可能なソフトウェアツールGAAS [ref.10]は、微生物ゲノムのデータベースに対するショットガン配列のBLAST検索に基づいてAGSを推定している。しかしメタゲノムデータおよびリファレンスゲノムの規模が大きくかつ増大しており、BLASTの速度が比較的遅いためこの方法は計算上実用的ではない。さらに、微生物群集は、しばしば、培養されていないかまたはゲノムシーケンスされていない、高い割合の「新規な」微生物を含む。よく研究されたヒトの腸内微生物でさえ、現在の微生物リファレンスゲノムでは、種の豊富さ（species abundance ）で43％と、richnessで58％は捕獲できないと推定されている[ref.16]。 GAASが新規分類群からなるメタゲノムについてAGSを正確に推定できるかどうかは不明である。 RaesらはBLASTXを用いて35の必須シングルコピー遺伝子セットにアサインされたリードの密度に基づいて平均ゲノムサイズを推定することでこの問題に対処する提案をした[ref.9]。この方法は著しく速くコミュニティの組成への依存度も低いが、300塩基対（bp）未満のリードについては設計されておらず、またソフトウェアはリリースされなかった。現在のシーケンシング技術は、より長いリードを生成し始めているが、膨大な量の既存のショートリードシーケンシングデータが残っており、例えば、MetaHIT [ref1,17-19]とHuman Microbiome Project（HMP）[ref.20]プロジェクトでは、ヒト共生微生物から平均96bpの長さで1300億以上のメタゲノミックリードを生成した[ref,21]。 Raes法が近代的なショートリードライブラリのAGSを推定できるかどうかは明らかではない。その他[ref.22,23]も同様の方法を記述しているが、これらは広範に検証されておらず、ソフトウェアとして利用可能になっていない。

　この問題に対処するために、著者らは、ショットガンメタゲノミックデータからAGSを迅速かつ正確に推定するMicrobeCensusを開発し、このツールを1,352のヒト共生微生物サンプルに適用した。 Raesらにも同様のアプローチを採用したが、有意な方法論的改善がなされ、50bpのショートリードを使用してAGSを正確に推定することができた。著者らは、AGSが、ヒト体内の位置でヒト共生微生物が著しく異なり、機能的および分類学的に大きな違いがあることを発見した。例えば、便のメタゲノムのAGSは、2.5〜5.8メガベース（Mb）の範囲であり、Bacteroidesの量と代謝、生合成、および二成分系に関連する遺伝子は正に相関したが、Firmicutes（フィルミクテス門）および膜輸送に関連する遺伝子はより小さいAGSを有するメタゲノムに豊富に存在していた。さらに著者らは、AGSが比較メタゲノム解析の主要なバイアスであり、正規化が異なる存在量の遺伝子検出を改善することを確認した。これらの発見はMicrobeCensusの開発なしにはできなかった。

 　MicrobeCensusはメタゲノムのシーケンスデータから微生物群集の平均ゲノムサイズ（ average genome size (AGS)）を推定するためのツール。AGSは、ほぼすべての細胞性微生物（バクテリア、アーキア、真菌）に存在する一組の普遍的な単一コピー遺伝子ファミリーにシーケンシングリードをアサインすることよって推定する。 これらの遺伝子はゲノムごとに1回だけ出現するため、微生物群集の平均ゲノムサイズは、これらの遺伝子にヒットしたリードの割合と反比例する。AGSが得られれば、サンプル中に存在する微生物ゲノムの総カバレッジを得ることが可能となる（genome equivalents = total bp sequenced / AGS in bp）。これは遺伝子の存在量を正規化するのに有用である。

f:id:kazumaxneo:20180711115509p:plain

MicrobeCensusのAGS推定ワークフロー。

 

 

MicrobeCensusに関するツイート。

 

簡単に言えば、どのような全細胞性生物でも厳密に１コピーだけ存在している遺伝子のセット（AGS）を定義し、その遺伝子のリードカウントを全リードカウントと比較してゲノムサイズを推定する方法論です。１コピーしか存在しない遺伝子なので、AGSのリードカウントはゲノムサイズによっては変化しません。そのため、ゲノムサイズが大きければ"whole ゲノムリードカウント/ AGSリードカウント"、は大きくなり、ゲノムサイズが小さければ"whole ゲノムリードカウント/ AGSリードカウント"は小さくなります。この性質を利用すれば、全リード数とAGSのリードカウントの比率からゲノムサイズは推定できます。この方法を複数ゲノムが混じったデータセットに適用すれば、得られるゲノムサイズ値は、複数ゲノムサイズの平均値になります。

 インストール 

mac os10.12のAnaconda3.4.2でテストした。

依存

  • Supported platforms: Mac OSX, Unix/Linux; Windows not currently supported
  • Python version 2 or 3
  • Numpy
  • BioPython 

本体　Github

#Anaconda環境ならcondaで導入できる（*1）。
conda install -c bioconda microbecensus

#Anacondaを使ってなくてもpipで導入できる。
pip install microbecensus 

*1 #condaでインストール時にpymummerとコンフリクトするとの警告があったのでpymummerをconda uninstallで消してから導入した。

> run_microbe_census.py -h

$ run_microbe_census.py -h

usage: run_microbe_census.py [-h] [-v] [-V] [-r RAPSEARCH] [-n NREADS]

                             [-t THREADS] [-e]

                             [-l {50,60,70,80,90,100,110,120,130,140,150,175,200,225,250,300,350,400,450,500}]

                             [-q MIN_QUALITY] [-m MEAN_QUALITY] [-d]

                             [-u MAX_UNKNOWN]

                             SEQFILES OUTFILE

 

Estimate average genome size from metagenomic data.

 

positional arguments:

  SEQFILES              path to input metagenome(s); for paired-end

                        metagenomes use commas to specify each file (ex:

                        read_1.fq.gz,read_2.fq.gz); can be FASTQ/FASTA; can be

                        gzip (.gz) or bzip (.bz2) compressed

  OUTFILE               path to output file containing results

 

optional arguments:

  -h, --help            show this help message and exit

  -v                    print program's progress to stdout (default = False)

  -V, --version         show program's version number and exit

  -r RAPSEARCH          path to external RAPsearch2 v2.15 binary; useful if

                        precompiled RAPsearch2 v2.15 binary included with

                        MicrobeCensus does not work on your system

 

Pipeline throughput (optional):

  -n NREADS             number of reads to sample from SEQFILES and use for

                        average genome size estimation. to use all reads set

                        to 100000000. (default = 2000000)

  -t THREADS            number of threads to use for database search (default

                        = 1)

  -e                    quit after average genome size is obtained and do not

                        estimate the number of genome equivalents in SEQFILES.

                        useful in combination with -n for quick tests (default

                        = False)

 

Quality control (optional):

  -l {50,60,70,80,90,100,110,120,130,140,150,175,200,225,250,300,350,400,450,500}

                        all reads trimmed to this length; reads shorter than

                        this discarded (default = median read length)

  -q MIN_QUALITY        minimum base-level PHRED quality score (default = -5;

                        no filtering)

  -m MEAN_QUALITY       minimum read-level PHRED quality score (default = -5;

                        no filtering)

  -d                    filter duplicate reads (default = False)

  -u MAX_UNKNOWN        max percent of unknown bases per read (default = 100

                        percent; no filtering)

——

 

テストラン

git clone https://github.com/snayfach/MicrobeCensus.git
cd /path/to/MicrobeCensus/test
python test_microbe_census.py

an 3 tests in 44.751s

 

OK

——

 

ラン

シーケンスデータはFASTQまたはFASTAフォーマット、gzip (.gz) か bzip (.bz2) 圧縮していても受け付ける。ペアエンドはコンマで区切って記載する。

run_microbe_census.py -t 12 -n 2000000 pair1.fq,pair2.fq output
  • -t     number of threads to use for database search (default= 1)
  • -n    number of reads to sample from SEQFILES and use for average genome size estimation. to use all reads set to 100000000. (default = 2000000)
  • -e    quit after average genome size is obtained and do not estimate the number of genome equivalents in SEQFILES. useful in combination with -n for quick tests (default = False)

コンタミネーションがかなりあったバクテリアのシーケンスデータの出力。

reads_sampled:  200000

trimmed_length: 300

min_quality:    -5

mean_quality:   -5

filter_dups:    False

max_unknown:    100

 

Results

average_genome_size:    7452243.21927

total_bases:    246786445

genome_equivalents:     33.1157260624

output (END)

 

importして他のメタゲノム解析パイプラインのモジュールとして使うこともできます。詳細はGithubで確認してください。

引用

Average genome size estimation improves comparative metagenomics and sheds light on the functional ecology of the human microbiome

Nayfach S, Pollard KS.

Genome Biol. 2015 Mar 25;16:51.