bcftoolsでvcfファイルの編集(個人メモ・追記予定)
※2023年7月19日 一部修正
VCF形式
詳しい説明は省略。GATK公式による解説がおそらく一番分かりやすい。 VCFファイルを整形したり、ほしい多型情報を抽出するツールとしてはbcftoolsとVCFtoolsが定番となっている。個人的にはbcftoolの方が、
- biocondaで配布している
- bamファイルからの多型抽出コマンドも備えていて、多型抽出からフィルタリングをシームレスに行える
という点で扱いやすいと感じた。
bcftoolsのインストール方法
anaconda環境があればbiocondaからインストールできる。
conda install -c bioconda bcftools
もしくはgithubからインストールする。
VCFの連結
bcftools concat を使う。引数-aをつけると最後のファイルので1行目を上書き。サンプル名、フィールド等が同じであればこれをつける。
# make list echo -e "file-1.vcf.gz\nfile-2.vcf.gz\nfile-3.vcf.gz" > vcflist # concatenate tables bcftools concat -a -f vcflist -o vcf_conc.vcf.gz # add index bcftools index -f vcf_conc.vcf.gz
VCFの情報を確認
# list samples bcftools query -l file.bcf # number of samples bcftools query -l file.bcf | wc -l # list of positions bcftools query -f '%POS\n' file.bcf
VCFのフィルタリング
-i '条件'を引数に取ることで、条件に合う多型サイトを抽出できる。条件に合うものを除外するには-e '条件'とする。 複数サンプルのmulti-vcfファイルの場合、一つでも条件を満たすサンプルがあれば、そのサイトは抽出される。全サンプルが満たす条件で抽出する場合には、-eを使う。 パイプで繋げる時は-Ou (output uncompressed BCF)をつけると早い。
(例) 1. FILTERフィールドがPASSであり(GTAKなどで事前にフィルターをかけていた場合)、かつQAUL (quality)が100以上であり、 2. 1~4番目の個体のうちgenotypeが0/0である個体が3個体以上の多型を抽出
bcftools view -O u -i 'FILTER="PASS" & QUAL>100' file.vcf.gz | bcftools view -i 'COUNT(FMT/GT[0-3]="0/0")>2' -o file_filtered.vcf
※フィールドのINFOとFMTの違い INFOは全サンプルを総合した数値。FMT(FORMAT)は個別サンプルにつけられる値を表す。たとえばINFO/DPは全サンプルのdepth合計だが、FMT/DPとすると1サンプルごとのdepthを示す。 ※&と&&の違い &で条件を繋げると、全条件が一つのサンプル内で満たされなければいけない。&&で繋げた場合、各条件がいずれかのサンプルで満たされていればよい。
VCFをテーブルに変換
bcftools query で必要な情報を取り出してテーブルに保存できる。-f '条件' の形式。-H:ヘッダーをつけて保存
(例) タブ区切りでCHROM, POS, REF, ALT, QUAL, AC, AFの値、その右にGT, AD, DPの値を個体ごとにタブで区切って出力、テキスト形式で保存。
bcftools query -H -f '%CHROM\t%POS\t%REF\t%ALT\t%QUAL\t%AC\t%AF[\t%GT %AD %DP]\n' file_filered.vcf.gz > file_filtered_GTADDP.txt
ジェノタイプ変換
setGTプライグインが提供されている。ただしデータの変換はvcfにあった情報を破壊することになるので慎重に行うべき。一度保存してから外部のツールを使うべきかもしれない。 (例) GQ (Genotype Quality: 個体ごとのGTの確からしさ)が20未満の場合、そのGTを./.に変換
bcftools +setGT input.vcf.gz -o transformed.vcf.gz -- -t q -i 'GQ<20' -n "./."
参考サイト
- 公式マニュアル。これを読み込めば十分。
- 有志によるチートシート
- +setGT pluginについて
http://samtools.github.io/bcftools/howtos/filtering.html
。これに関しては情報が少ない。