如何獲得比較可靠的變異位點？

拿到 VCF 文件以後，怎麼過濾才能得到比較可靠的變異位點呢？本文主要對 BWA + GATK 流程得到的 VCF 文件進行探討，且僅涉及 SNPs 位點的過濾、不涉及 INDELs。

VCF 文件的格式

VCF 文件的格式請看這裡。

VCF 文件的過濾

VCF 文件的過濾可以考慮以下幾個方面：

1) 過濾 INDEL 附近的 SNPs；

2) 提取 SNP 位點；

3) GATK hard filtering；

4) 過濾 DP、GQ；

5) 過濾 maf、missing；

6) 僅保留具有2個等位基因的位點；

7) 過濾掉高雜和度的位點。

這些過濾信息僅供大家參考，每個工作的側重點不同，所做的過濾也不一樣。接下來看一下每一步過濾如何操作：

1 過濾 INDEL 附近的 SNPs

即使運行了 Indel realignment，並且使用 haplotypecaller，GATK 也不能完美地處理 INDELs，換言之，INDELs 附近的 SNPs 可能是不準確的，可以考慮在分析之前過濾掉。

上述命令中 raw.vcf 為初始的未過濾的 VCF 文件；--SnpGap 代表過濾掉 INDEL 附近 5 bp 以內的 SNPs；-O v 表示輸出格式為 vcf；如果文件較大可以考慮設為 -O z，輸出壓縮的 VCF 文件； --output 用來指明輸出文件的名稱。

2 提取 SNP 位點

上面用兩種方式實現提取 SNP 位點的操作；

如果使用 vcftools，記得加上--recode-INFO-all 參數，這樣才能保留 INFO 的信息，後面還需要根據這列信息進行過濾。如果不加，整列信息將變為 . ！

3 GATK hard filtering

對於 GATK hard filtering 的參數，請參考官網的介紹；

這步運行後，GATK 並沒有直接刪掉不符合閾值的位點，僅僅是在 FILTER 那列標註 PASS 或標註每個過濾參數自定義的 filterName ；

在最初嘗試過濾參數組合的時候，給每個過濾參數的命名可參考：特有 (FS) + 共有 (_snp_filter) 的方式，以便查看每個過濾參數對數據的影響。

4 過濾 DP、GQ

GQ 代表 genotype quality，和鹼基質量值的表示一樣 ，20 表示錯誤率為 1%；DP 是過濾後的覆蓋度，一般 2 倍體不能低於 5， 單倍體不能低於 3， 如果測序深度很高，可以考慮用 10；需要注意的是：這步過濾僅僅是把不符合要求的個體基因型改為 ./. ，並沒有刪掉任何的 SNPs！這也是需要在過濾 maf、missing 之前過濾 DP、GQ 的原因。5 過濾 maf、missing

注意 vcftools 的 --max-missing 1.0 代表沒有 missing，設為 0.8 代表 missing 率為 20%；6 僅保留具有2個等位基因的位點

含有多個等位基因的位點也可能是同源基因造成的，過濾後可提高 SNPs 準確性；很多分析不支持含有 2 個以上的等位基因的位點；7 過濾掉高雜和度的位點如果一個位點的雜和度過高，很可能是旁系同源基因或其它因素的影響，需要進行過濾以保證 SNPs 的準確性。可以考慮使用 STACKS 軟體中的 populations 命令來過濾，或者自己寫一個簡單的腳本來過濾。如果喜歡我們的文章，歡迎訂閱我們的公眾號。

喜歡這篇文章嗎？立刻分享出去讓更多人知道吧！