还在为微生物重测序变异检测发愁？samtools帮助你！

• 下载地址 http://www.htslib.org/download/，利用wget下载得到samtools-1.5.tar.bz2，bcftools-1.5.tar.bz2；

• 解压缩 ： tar -jxvf bcftools-1.5.tar.bz2

tar -jxvf samtools-1.5.tar.bz2

• 安装 （详细的安装信息可查看bcftools-1.5(samtools-1.5)文件夹下的INSTALL）：

cd bcftools-1.5 ( samtools-1.5)

./configure –prefix=/where/to/install #设置安装路径

Make

Make install

• 安装完成后 ，我们后面所用到的可执行程序都在上面设置的安装路径（/where/to/install）的bin文件夹下；

• 配置环境变量

  a.若要临时修改环境变量，可直接在终端输入下面一行命令：

  Export PATH=/where/to/install/bin:$PATH

  b.要永久修改环境变量可将下面第一行添加到~/.bash_profile(针对当前用户)或者/etc/profile(针对所有用户)文件的末尾，再执行第二行命令即可：

  Export PATH=$PATH: /where/to/install/bin；

  source ~.bash_profile 或者source /etc/profile

现在samtools 和bcftools已经安装好了，下面我们就进行calling SNPs和Indels了。

1. 输入文件

比对结果文件sample.sort..bam和参考基因组序列ref.fa（上一期已对其建立索引，此处不再建立索引）;

2. Variant Calling主要流程

Samtools index sample.sort.bam #对bam文件建立索引,默认生成文件为bam文件加.bai,此处生成sample.sort.bam.bai;

Samtools rmdup sample.sort.bam sample.dedup.bam #去除PCR等实验过程中产生的多余duplications;

Samtools mpileup –q 20 –d 8000 –ugo sample.bcf –f ref.fa sample.dedup.bam #variant calls

参数说明：

-q,  --min-MQ     mapQ最小值，低于这个值则跳过；

-d,  --max-depth   最大覆盖深度；

-g,  --BCF        生成bcf格式文件；

-o,  --output FILE  输出文件；

-f ,  --fasta-ref FILE参考序列(fasta format)的索引文件名；

其他常用参数：

-A  --count-orphans  不丢弃异常配对reads；

-l   --positions FILE  包含区域位点的位置列表文件或者bed区域文件；

-r   --region REG    在指定区域产生pipeup，和-l一起使用；

-I    --skip-indel     不检测INDEL;

-m  --min-ireads 候选INDEL的最小间隔的reads数；

-v   --VCF          生成VCF格式文件；

-u  --uncompressed 产生未压缩的vcf/bcf文件；

Bcftools call –vmO z –o sample.vcf.gz sample.bcf #将bcf转化为vcf文件，bcf为vcf的二进制格式

参数说明：

-v   --variants-only       仅输出变异位点，不加此参数，则输入所有位点信息；

-c   --consensus-caller    the original calling method;

-m  -- multiallelic-caller   alternative model for multiallelic and rare-variant calling，calling的方法，与-c不能同时使用；

-O   --output-type 指定输出文件类型， 'b'表示压缩后BCF文件; 'u'表示未压缩的BCF文件; 'z'表示压缩的VCF文件; 'v' 表示未压缩的VCF文件；

bcftools filter -O z -o sample.filtered.vcf.gz -i 'QUAL>20 && DP>5' sample.vcf.gz #过滤质量值低于20并且深度小于5的变异位点

参数说明：

-O  --output-type 指定输出文件类型， 'b'表示压缩后BCF文件; 'u'表示未压缩的BCF文件; 'z'表示压缩的VCF文件; 'v' 表示未压缩的VCF文件；

-o,  --output FILE  输出文件；

-i   --Include 筛选出满足exp条件的变异位点；

-e     --exclude 剔除满足exp条件的变异位点，与-i相反；

--threads    线程数；

Note: 上面的过滤条件中可以加入DP4参数，过滤出高质量的变异位点

Vcf 文件格式：由#开头的注释部分（图1）和变异主体部分（图2，图3），注释部分主要包括软件版本，命令行，参考序列信息以及主体部分参数的注释，主体部分由以tab隔开的8列组成，依次为参考序列名，variant的位置，variant的ID，参考序列的碱基，Phred格式(Phred_scaled)的质量值，过滤结果和第8列variant 的详细信息（图3），采用tag=value形式，tag间用“;”隔开；

图1 VCF文件注释部分