2015-03-05
まず、
MiSeqデータのMHC領域へのマッピング(5):bwa, samtools indexまでを同時に行うスクリプト - NGSデータ解析まとめ
のshell script (mapping_mhc.sh)を使って、DRR003760, DRR003761, DRR003762の3つのMiSeqデータをヒトchr.6にマッピングし、.sort.bamファイルを作成して下さい。
ここでは、samtools mpileupコマンドを使用して、複数のmapping結果から、多型サイト(referenceのゲノム配列と異なる塩基サイト、indel)を抽出します。結果はBCFという形式のファイルに出力されます。その後、bcftoolsおよびvcfutils.plを用いて、VCF形式の多型サイトのリストを作成します。
Shell scriptを使用した場合、結果は./mapping_results/DRR0037*/の中にあるので、以下のようにコマンドを入力・実行します。
samtools mpileup -g -f Homo_sapiens.GRCh38.dna_rm.chromosome.6.fa ./mapping_results/DRR003760/DRR003760.sort.bam ./mapping_results/DRR003761/DRR003761.sort.bam ./mapping_results/DRR003762/DRR003762.sort.bam > MHC_test4.bcf
次に、以下のコマンドを実行します。ここでは・・・(<-- あとで調べること!)
bcftools view -bvcg MHC_test4.bcf > MHC_test4.var.bcf
次に、以下のコマンドを実行し、VCF形式のoutput fileを作ります(SNPsのリスト)。
bcftools view MHC_test4.var.bcf | vcfutils.pl varFilter -> MHC_test4.var.vcf