vcftools安装及基础用法

vcftools是一种可以对VCF文件和BCF文件进行格式转换及过滤的工具，功能非常强大，而且运算速度也很快。

1.下载及安装

1.1 下载地址

http://vcftools.sourceforge.net/downloads.html

1.2 安装

进入压缩包目录，进行解压。

$ cd /gss1/home/sunwei/software
$ tar xvf vcftools_0.1.13.tar.gz
$ cd /gss1/home/fzhang/sunwei/vcftools_0.1.13/bin

检查安装是否成功

$ vcftools
VCFtools (v0.1.13)
© Adam Auton and Anthony Marcketta 2009
Process Variant Call Format files
For a list of options, please go to:
        https://vcftools.github.io/examples.html
Questions, comments, and suggestions should be emailed to:
        vcftools-help@lists.sourceforge.net

2. 基础用法

2.1 vcf文件加ID

add id.pl拷贝至vcftools/bin目录下
add id.pl是一个老师写的脚本，这里不好直接放上来，所以需要添加id的话，请大家再去查找其他的教程，这里只是我自己做个备份。

perl add_id.pl root.hic.vcf root.hic.id.vcf

2.2 分开indel和SNP

只输出indel

vcftools --vcf  root.hic.id.vcf --keep-only-indels --recode --recode-INFO-all --out root.hic.id.indel

只保留SNP

vcftools --vcf root.hic.id.vcf --remove-indels --recode --recode-INFO-all --out root.hic.id.snp

2.3 vcf文件过滤

vcftools --vcf root.hic.id.vcf --max-missing 0.8 --maf 0.05 --min-alleles 2 --max-alleles 2 --recode --recode-INFO-all --out  root.hic.id.int0.8maf0.05.allele2

max-missing:分型完整度
maf：第二等位基因频率
min-alleles:最小等位基因个数
max-alleles:最大等位基因个数

哈迪温伯格平衡 hwe filtering（通常在人类中使用）

vcftools --vcf root.hic.id.vcf --remove-indels --max-missing 0.8 --maf 0.05 --min-alleles 2 --max-alleles 2 --hwe 0.01 --recode --recode-INFO-all --out root.hic.id.snp.hwe0.01

2.4 文件格式转换

转换vcf格式为ped map格式

vcftools --vcf root.hic.id.vcf --plink --out root.hic

2.5 vcf文件拆分

准备一个samplelist文件，即需要的样本的ID

vcftools --vcf hic.sort.ref.vcf --recode --recode-INFO-all --keep samplelist.txt --out root

详细说明书见官网：

http://vcftools.sourceforge.net/

2.6 统计等位基因频率

第一条染色体上的等位基因频率，注意vcf文件中chr的大小写。
结果中第一列是染色体；第二列SNP位置；第三列是这个位置一共出现了几种碱基，这里是两个；第四列是该位置出现的碱基总数，这里一个样本贡献了两个碱基位点；后面是该位置出现的碱基对应的频率。

$ vcftools --vcf root.id.vcf --freq --chr Chr1 --out Chr1_analysis
$ head -n 5 Chr1_analysis.frq
CHROM   POS     N_ALLELES       N_CHR   {ALLELE:FREQ}
Chr1    5154    2       286     C:0.713287      T:0.286713
Chr1    5187    2       296     A:0.614865      G:0.385135
Chr1    5220    2       282     A:0.241135      T:0.758865
Chr1    50889   2       294     A:0.870748      G:0.129252

计算整体文件中的等位基因频率

$ vcftools --vcf root.id.vcf --freq --out allel_analysis
$ tail -n 5 allel_analysis.frq
scaffold995     3028    2       304     G:0.927632      A:0.0723684
scaffold995     3082    2       304     G:0.842105      A:0.157895
scaffold995     3168    2       298     T:0.946309      C:0.0536913
scaffold995     3185    2       292     C:0.89726       T:0.10274
scaffold995     3228    2       302     A:0.407285      C:0.592715

2.7 比较两个vcf文件的变异位点

结果文件中的第一列是染色体；第二列和第三列是文件1和文件2中SNP位置；第四列中B表示两个文件中都有这个碱基，如果是1则表示只有文件1中有这个碱基，如果是2同理。

$ vcftools --vcf root.id.vcf --diff fei.id.vcf --diff-site --out in1_v_in2
$ head -n 5 in1_v_in2.diff.sites_in_files
CHROM   POS1    POS2    IN_FILE REF1    REF2    ALT1    ALT2
Chr1    5154    5154    B       C       C       T       T
Chr1    5187    5187    B       A       A       G       G
Chr1    5220    5220    B       A       A       T       T
Chr1    50889   50889   B       A       A       G       G

以上内容为个人粗浅的理解，欢迎讨论，如有错误，也敬请指出。引用请注明出处。