VCFtools

一、简介

官网:https://vcftools.github.io/index.html

说明: http://vcftools.sourceforge.net/perl_examples.html

教材: http://vcftools.sourceforge.net/man_latest.html

  • Filter out specific variants
  • Compare files
  • Summarize variants
  • Convert to different file types
  • Validate and merge files
  • Create intersections and subsets of variants

二、案例

语法:

vcftools [ --vcf FILE | --gzvcf FILE | --bcf FILE] [ --out OUTPUT PREFIX ] [ FILTERING OPTIONS ] [ OUTPUT OPTIONS ]

常用例子:

1.提取压缩的vcf中MT染色体上所有的位点的频率

vcftools --gzvcf /bioinfo/data/1000G/rawData/20130502/20130502.v5a/ALL.chrMT.phase3_callmom.20130502.genotypes.vcf.gz --freq --chr MT --out chrMT_analysis

2.去掉vcf中inde,并形成新的vcf

vcftools --gzvcf /bioinfo/data/1000G/rawData/20130502/20130502.v5a/ALL.chrMT.phase3_callmom.20130502.genotypes.vcf.gz --remove-indels --recode --recode-INFO-all --chr MT --out chrMT_snps_all

3.提取vcf中某个样品的信息

vcftools --gzvcf /bioinfo/data/1000G/rawData/20130502/20130502.v5a/ALL.chrMT.phase3_callmom.20130502.genotypes.vcf.gz --indv NA18603 --recode --out new.vcf

如果不加–recode,会报错Keeping individuals in ‘keep’ list

提取某些样本的信息

vcftools --vcf rename_rs_samples_filtered.vcf --remove skip_samples.tsv --recode --out new.vcf

skip_samples.tsv为一行一个样本

4.提取某个区间的信息

vcftools --gzvcf /bioinfo/data/1000G/rawData/20130502/20130502.v5a/ALL.chrMT.phase3_callmom.20130502.genotypes.vcf.gz --indv NA18603 --recode --bed /home/qqin/GenoNGS_dev/data/WellWise/wellwise.target.amplicon.bed --out NA18603.vcf

-bed –exclude-positions

3.合并vcf

样本的合并:

vcf-merge NA14624-100-1_gatk-ug.vcf.gz NA14624-100-1_gatk-ug.vcf.gz >merge_vcf

几个样本的合并,每个样本占有一列 位点的合并(可以多个样本的合并)

vcf-concat A.vcf.gz B.vcf.gz C.vcf.gz | gzip -c > out.vcf.gz
vcf-concat samples_ALL.chrX.phase3_shapeit2_mvncall_integrated.20130502.genotypes.vcf.recode.vcf samples_ALL.chr11.phase3_shapeit2_mvncall_integrated_v5.20130502.genotypes.vcf.recode.vcf >total.vcf

4.depth的计算

计算每个位点每个样本的测序深度

vcftools --geno-depth --vcf total_samples_gatk-ug.vcf

计算每个位点深度之和

vcftools --site-depth --vcf total_samples_gatk-ug.vcf

5.提取某一样本,某一个位点的信息

vcftools --snp rs7513464 --snp rs1058164 --indv HG22628 --gzvcf ALL.chr22.phase3_shapeit2_mvncall_integrated_v5a.20130502.genotypes.vcf.gz --recode --stdout

5.所有位点,用Annovar注释,然后提取

三、报错

1.tabix不能建立index

tabix -l NA14624-100-1_gatk-ug.vcf 报错:Could not load .tbi index of NA14624-100-1_gatk-ug.vcf

tabix -p vcf NA14624-100-1_gatk-ug.vcf 报错: Not a BGZF file

解决办法:

gunzip dbsnp_138.hg19.vcf.gz 
bgzip dbsnp_138.hg19.vcf 
tabix -p vcf dbsnp_138.hg19.vcf.gz

参考资料: 。。。

个人公众号,比较懒,很少更新,可以在上面提问题,如果回复不及时,可发邮件给我: tiehan@sina.cn

Sam avatar
About Sam
专注生物信息 专注转化医学