vcf 格式文件詳解


Vcf文件格式是GATK鍾愛的表示遺傳變異的一種文件格式。

就拿GATK給出的vcf例子說明吧,下面這個文件只表示了一個完整vcf文件的前幾個SNP。

這里寫圖片描述

看上去確實有點復雜,那就把它分為兩部分看吧,第一部分把他歸為說明文件,就是每一列最前面有2個#符號的那些列所提到的就是為了解釋下面“正文”INFO列中可能要出現的一些tags和和FORMAT列中對基因型的表示。第二部分可以歸為下面的內容:

#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT NA12878
chr1 873762 . T G 5231.78 PASS AC=1;AF=0.50;AN=2;DP=315;Dels=0.00;HRun=2;HaplotypeScore=15.11;MQ=91.05;MQ0=15;QD=16.61;SB=-1533.02;VQSLOD=-1.5473 GT:AD:DP:GQ:PL 0/1:173,141:282:99:255,0,255
chr1 877664 rs3828047 A G 3931.66 PASS AC=2;AF=1.00;AN=2;DB;DP=105;Dels=0.00;HRun=1;HaplotypeScore=1.59;MQ=92.52;MQ0=4;QD=37.44;SB=-1152.13;VQSLOD= 0.1185 GT:AD:DP:GQ:PL 1/1:0,105:94:99:255,255,0
chr1 899282 rs28548431 C T 71.77 PASS AC=1;AF=0.50;AN=2;DB;DP=4;Dels=0.00;HRun=0;HaplotypeScore=0.00;MQ=99.00;MQ0=0;QD=17.94;SB=-46.55;VQSLOD=-1.9148 GT:AD:DP:GQ:PL 0/1:1,3:4:25.92:103,0,26
chr1 974165 rs9442391 T C 29.84 LowQual AC=1;AF=0.50;AN=2;DB;DP=18;Dels=0.00;HRun=1;HaplotypeScore=0.16;MQ=95.26;MQ0=0;QD=1.66;SB=-0.98 GT:AD:DP:GQ:PL 0/1:14,4:14:60.91:61,0,255

CHROM: 表示變異位點是在哪個contig 里call出來的,如果是人類全基因組的話那就是chr1…chr22,chrX,Y,M了。

POS: 變異位點相對於參考基因組所在的位置,如果是indel,就是第一個鹼基所在的位置。

ID: 如果call出來的SNP存在於dbsnp數據庫里,就會顯示相應的dbsnp里的rs編號。

REF和REF: 在這個變異位點處,參考基因組中所對應的鹼基和研究對象基因組中所對應的鹼基。

QUAL: 可以理解為所call出來的變異位點的質量值。Q=-10lgP,Q表示質量值;P表示這個位點發生錯誤的概率。因此,如果想把錯誤率從控制在90%以上,P的閾值就是1/10,那lg(1/10)=-1,Q=(-10)*(-1)=10。同理,當Q=20時,錯誤率就控制在了0.01。

FILTER: 理想情況下,QUAL這個值應該是用所有的錯誤模型算出來的,這個值就可以代表正確的變異位點了,但是事實是做不到的。因此,還需要對原始變異位點做進一步的過濾。無論你用什么方法對變異位點進行過濾,過濾完了之后,在FILTER一欄都會留下過濾記錄,如果是通過了過濾標准,那么這些通過標准的好的變異位點的FILTER一欄就會注釋一個PASS,如果沒有通過過濾,就會在FILTER這一欄提示除了PASS的其他信息。如果這一欄是一個“.”的話,就說明沒有進行過任何過濾。

到現在,我們就可以解釋上面的例子了:

chr1:873762是一個新發現的T/G變異,並且有很高的可信度(qual=5231.78)。
chr1:877664是一個已知的變異為A/G 的SNP位點,名字rs3828047,並且具有很高的可信度(qual=3931.66)。
chr1:899282是一個已知的變異為C/T的SNP位點,名字rs28548431,但可信度較低(qual=71.77)。
chr1:974165是一個已知的變異為T/C的SNP位點,名字rs9442391,但是這個位點的質量值很低,被標
成了“LowQual”,在后續分析中可以被過濾掉。

Vcf文件看起來很復雜,挺嚇人的樣子,但是里面大部分都是一些tags,而這些tags基本上都是在VASR中過濾用的,能夠理解每個tags的意思最好,如果實在不理解也就不用管了。其實最關鍵的信息也就是那么幾列:

chr1 873762 . T G [CLIPPED] GT:AD:DP:GQ:PL 0/1:173,141:282:99:255,0,255
chr1 877664 rs3828047 A G [CLIPPED] GT:AD:DP:GQ:PL 1/1:0,105:94:99:255,255,0
chr1 899282 rs28548431 C T [CLIPPED] GT:AD:DP:GQ:PL 0/1:1,3:4:25.92:103,0,26

其中最后面兩列是相對應的,每一個tag對應一個或者一組值,如:

chr1:873762,GT對應0/1;AD對應173,141;DP對應282;GQ對應99;PL對應255,0,255。

GT: 表示這個樣本的基因型,對於一個二倍體生物,GT值表示的是這個樣本在這個位點所攜帶的兩個等位基因。0表示跟REF一樣;1表示表示跟ALT一樣;2表示第二個ALT。當只有一個ALT 等位基因的時候,0/0表示純和且跟REF一致;0/1表示雜合,兩個allele一個是ALT一個是REF;1/1表示純和且都為ALT;

AD: 對應兩個以逗號隔開的值,這兩個值分別表示覆蓋到REF和ALT鹼基的reads數,相當於支持REF和支持ALT的測序深度。

DP: 覆蓋到這個位點的總的reads數量,相當於這個位點的深度(並不是多有的reads數量,而是大概一定質量值要求的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。


注意!

本站转载的文章为个人学习借鉴使用,本站对版权不负任何法律责任。如果侵犯了您的隐私权益,请联系我们删除。



 
粤ICP备14056181号  © 2014-2021 ITdaan.com