macでインフォマティクス

macでインフォマティクス

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

bamからのリードの抽出と他のゲノムアセンブリへのリアライメントを素早く実行する Bazam

2019 4/20  論文引用

2020 11/12 誤字修正、タイトル修正

 

 過去10年間にわたるハイスループットゲノムシーケンシングマシンの大規模な採用は、巨大な可能性を有する膨大な量のゲノムデータを生み出してきた。ゲノムデータは、座標 (coordinate) ソートされたBAMまたはCRAMフォーマットでアライメントされたリードとしてよく保存されたり交換される。多くのアプリケーション(アライメントのviewingやルーチンのvariant callingなど)で直接使用できるため、この形式は一般的である。しかしながら、アライメントされた形式での保存は、データがリファレンスゲノムおよび使用されるアライメント方法に結び付けられるという重大な欠点を有する。多くの結果はこれらのパラメータに非常に敏感であり、組み合わされたデータセットは通常これらのパラメータが同一で​​ない限り全く一緒に分析することはできない。その結果、データを最適に利用するために、ユーザは最近ビルドされたゲノムおよびリファレンスにデータをリアライメントする必要があることが多い。これは、ゲノムデータを効率的にリアライメントするための能力に対する広範なそして増大する必要性をもたらしている。
 しかし、アライメントされたデータからペアになったリードを再位置合わせすることは、計算量が多く、標準的な方法では不便である。これらの問題は、アライナーが最適にアラインするためにペアリードの両方に同時にアクセスしなければならないために発生する。両方のリードは通常アライメントファイルに格納されているが、coordinateソートファイルでは、かなりの部分が互いに離れている可能性がある。これらの場合、ペアになった両方のリードを一緒に出力に書き込むことができるような高価なランダムルックアップが必要になる。したがって、標準的なリアライメント方法では、最初にすべてのリードを抽出し、次にリアライメント前にそれらをリード名でソートする。これは抽出を可能にするが、その過程は時間がかかりそして中間工程のためにかなりのコンピュータリソースを必要とする。興味深いことに、Picard Tools [ref.1]は別の方法、SamToFastqの形式でペアのリードを抽出する。これはペアのリードを抽出する際のこれら中間ステップの必要性を回避している。ただし、この方法はコミュニティでは広く使用されていない。  SamToFastqはメモリ使用のために最適化されておらず、大きなデータセットでの使用には実用的ではないためである。さらに、Picard Toolsは特定のlocusをターゲットにすることはできず、シングルの出力ストリームしか発行できないため、下流のシングルプロセスの最大スループット(アライメントなど)がボトルネックになる。
 ここでは、並列処理やその他の追加機能を提供しながら、メモリ使用を最適化するSamToFastqの代替手段であるBazamを紹介する。 Bazamは出力ストリームを別々のリアライメントのために複数のパスに分割することで並列性を高めている。この手法を使用すると、無制限の数の並列アライナを使用してソースのシングルアラインメントをリアライメントでき、計算クラスタまたはクラウドコンピューティングリソースが使用可能になったときのプロセスが大幅に高速化される。

 リアライメントは重要なアプリケーションだが、それ以外にBazamはペアリードの詳細な情報に依存する他のアプリケーションのためのユーティリティも提供する。応用例としては、クオリティ管理やstructural variant callingがある。 Bazamは、特に興味のある2つの追加機能を提供している:read position taggingとlocalised extractionである。Read position taggingは、ストリーミングされる各リード名を元のアライメントポジションが含まれるようにリネームする。この機能により、リアライメント後に新旧のアライメントアライメントポジションを簡単に比較できる。Localised extractionは、リアライメントが特定のゲノム座標と重複するリードに限定されることを可能にする。リアライメントと同様に、これは標準的なツールを使用して実現できる。ただし、これらのツールは、ペアの片方だけが関心領域に重なっている場合、ペアの両方のリードを出力することはしないため、ペアの両方のリードを必要とするアプリケーションには適していない。このPreprintではBazamの実装について説明し、それが精度を犠牲にすることなく効率を向上させることを実証する。

 

https://www.biorxiv.org/content/biorxiv/early/2018/10/04/433003.full.pdf

f:id:kazumaxneo:20190114161609p:plain

Different configurations for using Bazam.  Preprintより転載

 


 

 

インストール

ビルドするかjavaの実行ファイルをダウンロードする。またはcondaで導入する。

本体 Github

#Anaconda環境でcondaを使い導入
conda install -c bioconda -y bazam

> bazam -h

$ bazam -h

================================================================================

 

Bazam

 

================================================================================

 

error: Missing required option: bam

usage: java -jar bazam.jar -bam <bam> -L <regions>

 -bam <arg>           BAM file to extract read pairs from

 -dr <arg>            Specify a read name to debug: processing of the read

                      will be verbosey printed

 -f,--filter <arg>    Filter using specified groovy expression

 -gene <arg>          Extract region of given gene

 -h,--help            Show help

 -L,--regions <arg>   Regions to include reads (and mates of reads) from

 -n <arg>             Concurrency parameter (4)

 -namepos             Add original position to the read names

 -o <arg>             Output file

 -pad <arg>           Amount to pad regions by (0)

 -r1 <arg>            Output for R1 if extracting FASTQ in separate files

 -r2 <arg>            Output for R2 if extracting FASTQ in separate files

 -s <arg>             Sharding factor: format <n>,<N>: output only reads

                      belonging to shard n of N

 

This tool is built with Groovy NGS - the Groovy way to work with NGS data. 

 

実行方法

bwaを作って新しいリファレンスゲノムにマッピングし直す。 新しいリファレンスのindexingはあらかじめ行っておく。

bazam -bam my.bam | bwa mem -t 8 -p ref.fa - | \
samtools view -bSu - | \
samtools sort -o output.bam


#または
bazam -bam my.bam | bwa mem -t 8 -p ref.fa - | \
samtools sort -O BAM -@ 8 -o output.bam -
  • -bam   BAM file to extract read pairs from
  • -n    Concurrency parameter (4)

 

bazamだけ実行した場合はfastqが出力される。

bazam -bam test.bam > tmp.fastq

 

 特定の領域のアライメントだけ新しいリファレンスゲノムにマッピングし直す。 

bazam -bam my.bam -L chr1:5000000-6000000 | bwa mem -p ref.fa - | \
samtools view -bSu - | \
samtools sort -o out.bam
  • -L    Regions to include reads (and mates of reads) from.
     

 特定の領域のアライメントだけ新しいリファレンスゲノムにマッピングし直す。 4リードごとに1リードだけ(2リード目)だけ取り出す。

bazam -s 2,4 -bam my.bam | bwa mem -p ref.fa - | \
samtools view -bSu - | \
samtools sort -o out.bam
  • -L    Regions to include reads (and mates of reads) from.
     

 

Groovyの正規表現を使い、より高度なフィルタリングを行いながらfastqを取り出すことも可能になっている。詳細はGithubで確認して下さい。

引用

Bazam: a rapid method for read extraction and realignment of high-throughput sequencing data

Simon P. Sadedin, Alicia Oshlack
Genome Biology 2019 20:78

 

Bazam : A rapid method for read extraction and realignment of high throughput sequencing data

Simon P Sadedin, Alicia Oshlack

bioRxiv preprint first posted online Oct. 3, 2018

doi: https://doi.org/10.1101/433003

 

関連


 *1

ピークメモリ他の要素については調べてませんが、スピードだけで言えば、samtools fastq -@ threads の方が早く終わるかもしれません(本ツールは"130G or so of data could require up to 16G of RAM"という表現があるように、メモリ使用量が極めて少ない)。