VCF文件全称为Variant Call Format，通常为GATK和Samtools软件处理所得到。在没接触过的人乍一看确实很吓人，但是只要理解了其中的一些关键tag再进行过滤就可以有的放矢了，下面以如下例子做讲解。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
##fileformat=VCFv4.0
##fileDate=20090805
##source=myImputationProgramV3.1
##reference=1000GenomesPilot-NCBI36
##phasing=partial
##INFO=<ID=NS,Number=1,Type=Integer,Description="Number of Samples With Data">
##INFO=<ID=DP,Number=1,Type=Integer,Description="Total Depth">
##INFO=<ID=AF,Number=.,Type=Float,Description="Allele Frequency">
##INFO=<ID=AA,Number=1,Type=String,Description="Ancestral Allele">
##INFO=<ID=DB,Number=0,Type=Flag,Description="dbSNP membership, build 129">
##INFO=<ID=H2,Number=0,Type=Flag,Description="HapMap2 membership">
##FILTER=<ID=q10,Description="Quality below 10">
##FILTER=<ID=s50,Description="Less than 50% of samples have data">
##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">
##FORMAT=<ID=GQ,Number=1,Type=Integer,Description="Genotype Quality">
##FORMAT=<ID=DP,Number=1,Type=Integer,Description="Read Depth">
##FORMAT=<ID=HQ,Number=2,Type=Integer,Description="Haplotype Quality">
#CHROM POS     ID        REF ALT    QUAL FILTER INFO                              FORMAT      NA00001        NA00002        NA00003
20     14370   rs6054257 G      A       29   PASS   NS=3;DP=14;AF=0.5;DB;H2           GT:GQ:DP:HQ 0|0:48:1:51,51 1|0:48:8:51,51 1/1:43:5:.,.
20     17330   .         T      A       3    q10    NS=3;DP=11;AF=0.017               GT:GQ:DP:HQ 0|0:49:3:58,50 0|1:3:5:65,3   0/0:41:3
20     1110696 rs6040355 A      G,T     67   PASS   NS=2;DP=10;AF=0.333,0.667;AA=T;DB GT:GQ:DP:HQ 1|2:21:6:23,27 2|1:2:0:18,2   2/2:35:4
20     1230237 .         T      .       47   PASS   NS=3;DP=13;AA=T                   GT:GQ:DP:HQ 0|0:54:7:56,60 0|0:48:4:51,51 0/0:61:2
20     1234567 microsat1 GTCT   G,GTACT 50   PASS   NS=3;DP=9;AA=G                    GT:GQ:DP    0/1:35:4       0/2:17:2       1/1:40:3

VCF文件大致可分为三个部分：

  • 带有##的注释部分
  • 表头信息
  • INFO(变异详细信息)

注释部分

这一部分的信息做大致了解即可，通常包括VCF版本、参考基因组信息、运行命令以及介绍各个参数所代表的意义。此部分的格式分为四种：

1.key=value格式

1
2
3
4
5
##fileformat=VCFv4.0
##fileDate=20090805
##source=myImputationProgramV3.1
##reference=1000GenomesPilot-NCBI36
##phasing=partial

介绍了VCF版本，日期，参考基因组等信息
另一种为如下形式

2.INFO 

1
##INFO=<ID=ID,Number=number,Type=type,Description=”description”>

type通常为整数、浮点数、标识符、字符串

3.FILTERs

1
##FILTER=<ID=ID,Description=”description”>

4.FORMAT 

1
##FORMAT=<ID=ID,Number=number,Type=type,Description=”description”>

表头信息

表头信息分为以下：

  • CHROM

SNP或者INDEL所在染色体的名称

  • POS

SNP或者INDEL所在染色体的位置，通常INDEL为第一个碱基所在染色体的位置

  • ID

为SNP在dbSNP的编号，通常植物是没有的，显示为.

  • REF

此位置碱基在参考基因组的序列，通常为A、T、C、G、N里面的一种

  • ALT

此位置发生变异的碱基或者indel

  • QUAL

变异碱基的可信度，该值越高，变异的可信度越高;Q=-10lgP，Q表示质量值；P表示这个位点发生错误的概率。如果想把错误率从控制在90%以上，P的阈值就是0.1，那lg（0.1）=-1，Q=（-10）*（-1）=10。同理，当Q=20时，错误率就控制在了0.01。

  • FILTER

过滤情况，PASS说明质量较好，通过过滤；q10表示过滤质量低于10；s50表示低于总样本数的50%；'.'代表没有经过任何过滤

例子

20 14370 rs6054257 G A 29 PASS
表示这个变异位于20号染色体，14370这个碱基位置，在SNP数据库中对应ID为rs6054257，参考基因组碱基为G，变异碱基为A，质量值为29，通过过滤

20 17330 . T A 3 q10
表示这个变异位于20号染色体，17330这个碱基位置，在SNP数据库中无对应ID，参考基因组碱基为T，变异碱基为A，质量值为3，未通过过滤

20 1110696 rs6040355 A G,T 67 PASS
表示这个变异位于20号染色体，1110696这个碱基位置，在SNP数据库中对应ID为1110696，参考基因组碱基为A，变异碱基为A或T，质量值为67，通过过滤

变异详细信息

  • GT 表示这个样本的基因型，对于一个二倍体生物，GT值表示的是这个样本在这个位点所携带的两个等位基因。0表示跟REF一样；1表示表示跟ALT一样；2表示第二个ALT。
    如20 14370 rs6054257 G A 29 PASS
    0/0–>G/G
    0/1–>G/A
    1/1–>A/A
  • AD 这两个值分别表示覆盖到REF和ALT碱基的reads数，相当于支持REF和支持ALT的测序深度。
  • DP 覆盖到这个位点的总的reads数量，相当于这个位点的深度.
  • PL 对应3个以逗号隔开的值，这三个值分别表示该位点基因型是0/0，0/1，1/1的没经过先验的标准化Phred-scaled似然值（L）。如果转换成支持该基因型概率（P）的话，由于L=-10lgP，那么P=10^（-L/10），因此，当L值为0时，P=10^0=1。因此，这个值越小，支持概率就越大，也就是说是这个基因型的可能性越大。
  • GQ 表示最可能的基因型的质量值。表示的意义同QUAL。

chr1 899282 rs28548431 C T [CLIPPED] GT:AD:DP:GQ:PL 0/1:1,3:4:25.92:103,0,26
在这个位点，GT=0/1，也就是说这个位点的基因型是C/T；GQ=25.92，质量值并不算太高，可能是因为cover到这个位点的reads数太少，DP=4，也就是说只有4条reads支持这个地方的变异；AD=1,3，也就是说支持REF的read有一条，支持ALT的有3条；在PL里，这个位点基因型的不确定性就表现的更突出了，0/1的PL值为0，虽然支持0/1的概率很高；但是1/1的PL值只有26，也就是说还有10^(-2.6)=0.25%的可能性是1/1；但几乎不可能是0/0，因为支持0/0的概率只有10^(-10.3)=5*10-11。

  • AC AC(Allele Count) 表示该Allele的数目
  • AF(Allele Frequency) 表示Allele的频率
  • AN(Allele Number) 表示Allele的总数目

对于1个diploid sample而言：则基因型 0/1 表示sample为杂合子，Allele数为1(双倍体的sample在该位点只有1个等位基因发生了突变)，Allele的频率为0.5(双倍体的 sample在该位点只有50%的等位基因发生了突变)，总的Allele为2； 基因型 1/1 则表示sample为纯合的，Allele数为2，Allele的频率为1，总的Allele为2。

  • Dels Fraction of Reads Containing Spanning Deletions。进行SNP和INDEL calling的结果中，有该TAG并且值为0表示该位点为SNP，没有则为INDEL。
  • FS 使用Fisher’s精确检验来检测strand bias而得到的Fhred格式的p值。该值越小越好。一般进行filter的时候，可以设置 FS < 10～20。
  • MQ 表示覆盖序列质量的均方值RMS Mapping Quality
  • FQ phred值关于所有样本相似的可能性
  • AF1 AF(Allele Frequency) 表示Allele的频率，AF1为第一个ALT allele 发生频率的可能性评估
  • AC1 AC表示Allele（等位基因）的数目,AC1为对第一个ALT allele count的最大可能性评估
  • IS 插入缺失或部分插入缺失的reads允许的最大数量
  • G3 ML 评估基因型出现的频率
  • HWE chi^2基于HWE的测试p值和G3
  • CLR 在受到或者不受限制的情况下基因型出现可能性log值
  • UGT 最可能不受限制的三种基因型结构
  • CGT 最可能受限制三种基因型的结构
  • PV4 四种P值得误差，分别是（strand、baseQ、mapQ、tail distance bias）
  • INDEL 表示该位置的变异是插入缺失
  • PC2 非参考等位基因的phred（变异的可能性）值在两个分组中大小不同
  • PCHI2 后加权chi^2，根据p值来测试两组样本之间的联系
  • QCHI2 Phred scaled PCHI2.
  • PR 置换产生的一个较小的PCHI2
  • QBD Quality by Depth，测序深度对质量的影响
  • RPB 序列的误差位置（Read Position Bias）
  • MDV样本中高质量非参考序列的最大数目
  • VDBVariant Distance Bias，RNA序列中过滤人工拼接序列的变异误差范围

参考资料
http://www.internationalgenome.org/wiki/Analysis/vcf4.0/
http://samtools.github.io/hts-specs/VCFv4.2.pdf
http://samtools.github.io/bcftools/bcftools.html
http://www.cnblogs.com/emanlee/p/4562064.html
http://www.bio-info-trainee.com/863.html