vcfファイルのサンプル名がbamファイル名なんですけど…
bcftoolsで何も考えずにバリアントコールをしたとき、vcfのサンプル名はbamファイル名になっています。 別にこのままでもいいかもしれませんがなんかダサい感じがします。
そこでvcfの列名の部分にサンプル名が来るようにしたいと思います。バリアントコールは以下のように行う ことにします。
bcftools mpileup -f reference.fa sample1.bam sample2.bam sample3.bam | bcftools call -mv -Ob -o calls.bcf
ここで、vcfファイルのヘッダーに対し、任意のサンプル名をつけるためには bcftools mpileup
コマンドにおいて追加で -G
または --read-groups
というオプションを用います。
-G
オプションはファイル名を引数にとり、このファイルにbamファイルとサンプル名の対応関係を記録しておくこと
でその情報を bcftools mpileup
にわたします。
この際に指定するファイルの形式は以下のようになります(ただし、1つのbamファイルに1つのサンプル由来の リードのみ含まれているとします)。
* sample1.bam sample1 * sample2.bam sample2 * sample3.bam sample3
1列目には*
を書きます。これは1つのbamファイルに含まれるリードはすべて1つのサンプル由来であることを
表しています。もしそうでない(1つのbamファイルに複数サンプルのリードが含まれる)ときは
リードグループIDを1列目に指定することになります(詳しくは公式マニュアル参照)。
続いて2列目にはbamファイルのファイル名を指定し、
それぞれに対応するサンプル名を3列目に記載するといったファイル形式です。区切り文字は空白を用います。
上記のファイルを read_groups.txt
としてカレントディレクトリに保存した場合、具体的な
コマンドラインは以下のようになります。
bcftools mpileup -f reference.fa -G read_groups.txt sample1.bam sample2.bam sample3.bam | bcftools call -mv -Ob -o calls.bcf
これでvcfファイルの列名がちゃんとサンプル名になったと思います。
すでにあるvcfのサンプル名を変えたいんですけど…
もしかしたら、もうすでにあるvcfのサンプル名を変えたいといった状況に遭遇するかもしれません。
その場合は bcftools reheader
が有用です。-s
または --samples
オプションを用います。
-s
オプションも同様にファイル名を引数にとります。このファイルには各行ごとにサンプル名を
記載するだけで大丈夫です。このとき、vcfの列順とファイルに記載するサンプル名の順が一致していること
が必要になります。
以下に例を示します。
sample1 sample2 sample3
上記のようなファイルを samples.txt
というファイル名でカレントディレクトリに保存したら、
以下のコマンドを実行します。
bcftools reheader -s samples.txt -o new.bcf old.vcf
次回はこのあたりを踏まえ、マッピングからバリアントコールまでのパイプラインを snakemake
を用いて
記述していきたいと思います。