通用過濾
Vcftools(http://vcftools.sourceforge.net) 對vcf文件進行過濾
第一步:過濾最低質量低於30,次等位基因深度(minor allele count)不少於3
1 vcftools --gzvcf raw.vcf.gz --max-missing 0.5 --mac 3 --minQ 30 --recode --recode-INFO-all --out raw.g5mac3 2 3 ## 4 --gzvcf 壓縮的vcf文件; 5 --max-missing 過濾掉缺失率大於50%的位點; 6 --minQ 過濾掉低於30的質量粉豬; 7 --Mac 次要等位基因深度為3,過濾小於3的位點; 8 --recode 輸出過濾后的VCF文件 9 --recode-INFO-all 包含原來文件中所有的INFO信息 10 --out 輸出文件
第二步:上述結果文件raw.g5mac3.recode.vcf, 基於最低深度進行過濾
vcftools --vcf raw.g5mac3.recode.vcf --minDP 3 --recode --recode-INFO-all --out raw.g5mac3dp3 ## --minDP 最低的深度
第三步:刪除缺失率過高的樣本
1 ## 查看各個樣本的缺失率 2 vcftools --vcf raw.g5mac3dp3.recode.vcf --missing-indv 3 --missing-indv 查看每個樣本的variant缺失情況。輸出文件"out.imiss",最后一列即表示樣本中的variants缺失率; 4 5 6 ## 利用awk 輸出缺失率》0.5的樣本 7 awk '($5 >0.5){print $0}' out.imiss |cut -f1 >lowDP.indv 8 9 ## 去除高缺失率的樣本 10 vcftools --vcf raw.g5mac3dp3.recode.vcf --remove lowDP.indv --recode --recode-INFO-all --out raw.g5mac3dplm 11 12 --remove 根據lowDP.indv列表去除vcf文件中的高缺失率樣本
第四步:基於最大缺失比例,平均深度和次等位基因頻率(MAF)過濾
1 vcftools --vcf raw.g5mac3dplm.recode.vcf --max-missing 0.95 --maf 0.05 --recode --recode-INFO-all --out DP3g95maf05 --min-meanDP 20 2 3 ## 4 --max-missing 最大缺失率 5 --maf 次等位基因頻率,排在第二位的基因頻率
若你的群體來此多個區域,你想對不同的群體的樣本進行分布過濾,可進行如下操作
1 ## 提供一個樣本信息popmap 2 BR_002 BR 3 BR_004 BR 4 BR_006 BR 5 BR_009 BR 6 BR_013 BR 7 BR_015 BR 8 BR_016 BR 9 BR_021 WL 10 BR_023 WL 11 BR_024 WL 12 13 ## 根據第二列信息進行拆分 14 awk '($2=='BR'){print $0}' popman >1.keep 15 awk '($2=='WL'){print $0}' popman >2.keep 16 17 ## 用VCFtools分別估計不同群體的缺失比例 18 vcftools --vcf DP3g95maf05.recode.vcf --keep 1.keep --missing-site --out 1 19 vcftools --vcf DP3g95maf05.recode.vcf --keep 2.keep --missing-site --out 2 20 21 ## 合並兩個文本,根據最后一列提出缺失率大於0.1的樣本 22 cat 1.lmiss 2.lmiss | mawk '!/CHR/' | mawk '$6 > 0.1' | cut -f1,2 >> badloci 23 24 ## 利用VCFtools進行過濾 25 vcftools --vcf DP3g95maf05.recode.vcf --exclude-positions badloci --recode --recode-INFO-all --out DP3g95p5maf05
SNP過濾教程(一)