Bioinfo|bedtools-操作VCF文件

简介: Bioinfo|bedtools-操作VCF文件

   初步设想在Bioinfo板块中分享一些常见的生信分析软件的使用,原则就是有现成的轮子就不去自己造了。


   当我们得到原始的测序数据(raw data)后,经过初步的QC质控得到clean data (公司或者自己的),可以通过FastQC软件对结果进行一系列的统计和图形展示R|fastqcr QC数据处理。当确定数据得到要求后,进行mapping分析得到bam文件后,就是重头戏检测SNP了。

bcftools — utilities for variant calling and manipulating VCFs and BCFs 。


   呐,文本就通过几个常见的参数简单介绍一下bedtools是怎么“manipulating” VCF文件的。


一 mpileup 函数

使用方法:bcftools mpileup [OPTIONS] -f ref.fa in.bam [in2.bam […]]


示例:bcftools mpileup -C 50 -m 2 -F 0.002  -d 1000 -q 20 -Q 20 -a DP,DP4,ADR -f hg19.fa A.bam B.bam

参数说明:

-f: 指定参考基因组,后跟需要变异检测的bam文件,bam较少时可以直接写;

bcftools mpileup  -f hg19.fa A.bam B.bam

-b:后跟一个bamlist的文件;

bcftools mpileup -f hg19.fa -b bamList

-C --adjust-MQ INT  矫正的MQ值,推荐是50;

-q, -min-MQ INT (MQ质量值)# 过滤MQ质量值;

-Q, --min-BQ INT (base质量值)# 过滤base质量值;

-r, --regions 只call特定的染色体,或者点,区域

形式为CHR|CHR:POS|CHR:FROM-TO|CHR:FROM-[,…];


示例:
bcftools mpileup -r chr1 -f hg19.fa  A.bam | bcftools call -c -v > A.chr1.vcf
bcftools mpileup -r chr1:111111-122222 -f hg19.fa  A.bam B.bam | bcftools call -c -v > A.chr1_111111-122222.vcf

-R, --regions-file FILE  当有多个region时,将region存入文件,使用-R参数,文件格式:tab分割的三列,chr start end即可 ;

-a, --annotate LIST    在INFO/FORMAT列中增加一些如DP,DP4,AD等的信息,群call可能会需要;

Comma-separated list of FORMAT and INFO tags to output.


二 bcftools call函数

使用方法:bcftools call [OPTIONS] FILE


示例:bcftools mpileup  -f hg19.fa A.bam B.bam | bcftools call -c -v -r chr1:111111-122222 -s A -V indels

示例说明:

只检测A样本,且只输出chr1:111111-122222区域的SNP信息,保存到test2.vcf文件中。


参数说明:

-c  -m为两种call的方式,旧版本使用-c ,现在默认使用-m,能解决大部分任务

-V, --skip-variants snps|indels 跳过SNP或者INDEL位点

-v, --variants-only  只输出变异位点

-s, --samples LIST 只检测此处给出的样本ID (通用参数)

-S, --samples-file FILE 只对此文件中列出的样本进行检测 (通用参数)

-r  -R为通用参数,与上面用法一致


三 bcftools filter 函数

使用方法:bcftools filter [OPTIONS] FILE


示例:bcftools filter  ALL.vcf TYPE="snp" -e 'DP < 20 || MQ < 20 || %QUAL<10 || DP < 20 ||MAF < 0.05 || '

示例说明:群体检测结果ALL.vcf文件中,只保留 质量值大于10,MQ大于20,MAF大于等于0.05,样本总深度大于20X的过滤后的SNP数据集。


参数说明:

TYPE="snp":只保留SNP信息;

-s:样本过滤,bcftools filter ALL.vcf -s A

-e 可加感兴趣的各种参数,常见的如下:

群体样本总深度DP:bcftools filter  ALL.vcf TYPE="snp" -e 'DP < 20'

前两个样本的DP:bcftools filter  ALL.vcf -e 'FORMAT/DP[0-1] < 20 '

最小等位基因频率:bcftools filter  ALL.vcf -e 'MAF < 0.05 ' |ee

质量值和mapping 质量值:bcftools filter  ALL.vcf -e ‘%QUAL<10 || MQ < 20 '

BAF值(alt depth / (ref depth + alt depth ))


熟练使用以上及链接中的参数,就可以办到不写py或者pl脚本,强势过滤了。


四 bcftools index 函数

bgzip 压缩 vcf 文件为 gz 文件


bgzip -c A.vcf >A.vcf.gz ;bgzip -c B.vcf >B.vcf.gz

bcftools 为 gz 文件建索引

bcftools index -t A.vcf.gz ; bcftools index -t B.vcf.gz

五  bcftools merge 函数

使用说明:将多个VCF文件,合并成一个VCF文件

bcftools merge -m snps -f PASS,. --force-samples A.vcf.gz B.vcf.gz > A_B.merge.vcf

同样不用自己写脚本合并VCF文件,省事 且较少出错。


   建议多测试,熟悉以上参数的意义以及结果的不同,找到符合自己需求的参数组合,以得到理想的结果。


好了 ,拿着bedtools 去操纵VCF文件吧。

相关实践学习
RocketMQ一站式入门使用
从源码编译、部署broker、部署namesrv,使用java客户端首发消息等一站式入门RocketMQ。
消息队列 MNS 入门课程
1、消息队列MNS简介 本节课介绍消息队列的MNS的基础概念 2、消息队列MNS特性 本节课介绍消息队列的MNS的主要特性 3、MNS的最佳实践及场景应用 本节课介绍消息队列的MNS的最佳实践及场景应用案例 4、手把手系列:消息队列MNS实操讲 本节课介绍消息队列的MNS的实际操作演示 5、动手实验:基于MNS,0基础轻松构建 Web Client 本节课带您一起基于MNS,0基础轻松构建 Web Client
相关文章
|
4天前
|
弹性计算 运维 Shell
批量将txt 文件修改为doc 文件
【4月更文挑战第29天】
20 1
|
9月前
|
Web App开发 存储 数据可视化
VCF文件修改染色体格式的方法
VCF文件修改染色体格式的方法
|
10月前
ENVI_IDL:如何读取HDF4文件和HDF5文件的数据集和属性?
参考Documentation Center (l3harrisgeospatial.com) 由于时间和精力,这里并没有对其中的函数进行详细的说明,只是介绍了函数作用以及用法。
136 0
|
10月前
ENVI_IDL: 文本文件的读取(主要是txt、csv文件)
ENVI_IDL: 文本文件的读取(主要是txt、csv文件)
304 0
|
10月前
|
存储
ENVI_IDL: 创建HDF5文件并写入数据(以将Geotiff文件写入HDF文件为例) + 详细解析
ENVI_IDL: 创建HDF5文件并写入数据(以将Geotiff文件写入HDF文件为例) + 详细解析
120 0
|
10月前
|
存储 编解码
ENVI_IDL:读取OMI数据(HDF5)并输出为Geotiff文件+详细解析
ENVI_IDL:读取OMI数据(HDF5)并输出为Geotiff文件+详细解析
287 1
|
10月前
ENVI_IDL:(批量处理)如何对HDF5文件进行GLT文件的创建并进行几何校正最终输出为IMG格式?
ENVI_IDL:(批量处理)如何对HDF5文件进行GLT文件的创建并进行几何校正最终输出为IMG格式?
108 0
|
11月前
|
Python
python 读取 .csv/.xlsx/.xls 文件数据—批量修改文件名
python 读取 .csv/.xlsx/.xls 文件数据—批量修改文件名
|
11月前
|
数据处理
R|批量循环处理同一格式文件-csv,txt,excel
R|批量循环处理同一格式文件-csv,txt,excel
|
11月前
|
算法 开发者