利用SHAPEIT将vcf文件进行基因型(genotype)定相(phasing):查看两个突变是否来源于同一条链(染色体或父本或母本),two mutations carried by the same read

首先,下载SHAPEIT.

按照里面的步骤安装完后,将vcf文件进行基因型定相,分四步走。

第一步,将vcf文件转化为plink二进制文件(.bed, .bim, .fam)。

这一步需要用到GATK里的GenomeAnalysisTK工具,见如下命令:

java -Xmx8g -jar GenomeAnalysisTK.jar -T VariantsToBinaryPed  -R GRCh37.fa -V file.vcf --metaData sampleID.fam -mgq 0 -bed file.bed -bim file.bim  -fam file.fam 

 注:这里的metaData,输入的fam文件跟plink的fam格式一样,如果你没有fam文件,则需要自己手动生成,第一列和第一列可以都为样本的ID;

FAM files

The fields in a FAM file are

  • Family ID
  • Sample ID
  • Paternal ID
  • Maternal ID
  • Sex (1=male; 2=female; other=unknown)
  • Affection (0=unknown; 1=unaffected; 2=affected)
Example of a FAM file of the binary PLINK format:
FAM1 NA06985 0 0 1 1
FAM1 NA06991 0 0 1 1
0 NA06993 0 0 1 1
0 NA06994 0 0 1 1
0 NA07000 0 0 2 1
0 NA07019 0 0 1 1
0 NA07022 0 0 2 1
0 NA07029 0 0 1 1
FAM2 NA07056 0 0 0 2
FAM2 NA07345 0 0 1 1

######号外,在这一步的时候,如果出现报错信息:ERROR MESSAGE: No metadata provided for sample XXXXXX,很有可能是fam文件在作怪,需要对fam文件进行样本排序,使其的顺序跟vcf文件保持一致;

第二步,去除高缺失率的基因型和样本

plink --noweb --bfile file --keep-allele-order --me 1 1 --set-me-missing --make-bed --out file.nomendel
plink --noweb --bfile file.nomendel --keep-allele-order --geno 0.05 --make-bed --out file.nomendel.filter

  

 ####对数据进行高缺失率去除是很有必要的,如果跳过这一步,直接定相的话,会出现:“ERROR: XXX fully missing individuals (=100%)” 报错 

第三步,对基因型进行定位,需要下载遗传图谱(genetic map

shapeit --duohmm -W 5 --thread 1 --input-bed file.nomendel.filter.bed file.nomendel.filter.bim file.nomendel.filter.fam --input-map genetic_map_b37.txt -O hapData

第四步,将定相的文件生成vcf格式:

shapeit -convert --input-haps hapData --output-vcf hapData.vcf

  至此,vcf中没有定相的变异会由“0/0 0/1 1/0 1/1 ”变为“0|0 0|1 1|0 1|1”

  

原文地址:https://www.cnblogs.com/chenwenyan/p/9204535.html