アメリエフのブログ

Amelieff Staff Blog

bamにread groupを追記する

GATKは、BAMのフォーマットに厳しく(参照ページ)、たとえばヘッダにサンプル名を含むread groupのリストがあり、かつすべてのリードがそのread groupに属しているBAMしか受け付けません。

Read group(以下RG)は、たとえばBWAではマッピングのときに -R で指定することができます。

$ bwa mem -M -R "@RG¥tID:sample¥tSM:sample¥tPL:Illumina" genome.fa sample_R1.fastq.gz sample_R2.fastq.gz


マッピングの時にRGを付けていなくても、Picard toolkitのAddOrReplaceReadGroupsというツールで追記することができます。

inputのsample.bamにRGを追記してoutput.bamを作成する場合の例を見てみましょう。

$ java -jar AddOrReplaceReadGroups.jar I=sample.bam O=output.bam RGID=sample RGLB=sample RGPL=Illumina RGPU=run_barcode RGSM=sample


AddOrReplaceReadGroups.jarの必須オプションは以下の通りです。

  1. INPUT=File
  2. OUTPUT=File
  3. RGLB=String
  4. RGPL=String
  5. RGPU=String
  6. RGSM=String


Picardに比べるとだいぶ手間と時間がかかりますが、samtools mergeで追記する方法もあります(参照ページ)。
ご参考までに、手順は以下の通りです。

$ perl -e 'print "\@RG\tID:sample\tSM:sample\tLB:sample\tPL:Illumina\n"' > rg.txt
$ cat rg.txt
@RG     ID:sample       SM:sample       LB:sample       PL:Illumina
$ samtools merge -rh rg.txt output.bam sample.bam

【追記】 samtoolsにreheaderという機能が追加されたので、merge機能を使わなくても簡単にread groupを置換できるようになりました。