2018awk6回答例

1.5日目練習問題で使用したexample.vcfファイルの10列目、11列目にはそれぞれ「sample1」、「sample2」の検体の変異情報が記載されている。

「:」区切りの最初の1項目目がジェノタイプ情報であり、下記の意味になる。
0/0 : 変異なし
0/1 : 片方(父親由来もしくは母親由来のどちらか)の染色体に変異あり(ヘテロ)
1/1 : 両方の染色体に変異あり(ホモ)
./.    : デプスが不十分で判定不可能

ここで数字の1が2、3、4・・・などになることがあるが、それはVCFの5列目(ALT)の変異した配列の何番目になるかを示している。例:ALT→A,Cのとき、0/1であれば、Aの変異を1つ、0/2であればCの変異を1つ持つ。

sample1でデプスが不十分で判定不可能とされる変異はいくつか?sample2では?

#sample1
grep -v "^#" example.vcf|awk -F'\t' '{if($10~"^[.]/[.]"){cnt++}} END{print cnt}'

#sample2
grep -v "^#" example.vcf|awk -F'\t' '{if($11~"^[.]/[.]"){cnt++}} END{print cnt}'

また、DP(各検体のリード数)が4以下の変異数はsample1, 2それぞれ幾つか?

#sample1
grep -v "^#" example.vcf|awk -F'\t' '
 {
  n++;
  split($10,arr,":");
  if(arr[3]<=4){cnt++}
 }
 END{print n, cnt}
'

2.4日目で使用した、take.blastn.txtにはRで読み込むと意図しない結果になる文字(#)やほかのプログラムでエラーになりそうな文字が含まれている。次のそれぞれの文字が含まれている行数を出し、すべて「_」に置換せよ。

# ' ( )

awk -F'\t' '
{
 if($0~"[#'"'"'()]"){cnt++};
 gsub("[#'"'"'()]","_",$0);
 print $0
} 
END{print cnt}
' take.blastn.txt

3.とあるプログラムを使うためにtake.blastn.txtの2列目の遺伝子名を「“」で囲む必要が出てきた。2列目を「”」で囲んで出力せよ。

awk -F'\t' '
{
 OFS="\t"; #$2を変更するともともとタブ区切りだったのがスペースに変換されてしまうので、タブのまま出力するのに必要
 $2="\""$2"\"";
 print $0
}
' take.blastn.txt
  • 2018awk6回答例.txt
  • 最終更新: 2019/07/15 02:46
  • by 127.0.0.1