美文网首页生物信息学
10.18 bcftools下载使用 samtools语句

10.18 bcftools下载使用 samtools语句

作者: KK_f2d5 | 来源:发表于2018-10-18 16:26 被阅读0次

    师兄让下载bcftools,学习一下。

    1.下载

    电脑太不给力了,wget都没有,用brew下载一个

    太慢了,在github直接下载。

    https://github.com/samtools/bcftools/releases

    安装看readme

    共享目录出错,

    这个时候就要在/etc/ld.so.conf中加入xxx.so所在的目录。

    一般而言,有很多so档会在/usr/local/lib这个目录下

    在/etc/ld.so.conf中加入/usr/local/lib这一行后,执行 在etc中执行ldconfig 命令来更新生效。

    ./bftools 使用


    2.samtools语句

    介绍samtools和bcftools的语句与用法 :https://www.jianshu.com/p/da32e3c3a168

    samtools基本语句和功能:

    1.view 查看bam,sam文件的内容

    2.sort用来对bam文件进行排序

    3.merge将2个或2个以上的已经sort的bam文件合并成一个bam文件,合并后的文件不需要再次sort,是已经sort过了的。

    4.index对bam文件建立索引,生成后缀为.bai的文件,用于快速检索reads。需要对bam文件先进行排序,否则会报错。很多时候都需要索引文件的存在,特别是显示序列比对的情况下。比如samtools tview,gbrowse2等。

    5.faidx对fasta序列建立索引,生成后缀为.fai的文件。该命令也能依据索引文件快速提取fasta文件中的某一条序列。

    6.tview能直观的显示出reads比对到基因组的情况,和基因组浏览器有点类似。

    7.flagstat给出bam文件的比对结果的summary。

    8.depth得到每个位点的测序深度,并输出到标准输出

    samtools depth a.bam | less

    9.rmdup去除PCR重复和光学重复,重复的reads仅保留一对。

    10.reheader替换bam文件的头文件

    11.cat链接多个bam文件,适用于非sorted的bam文件

    12.mpileup(现在也属于bcftools)

    使用samtools的子命令mpileup分析参考序列上的每个碱基位点的比对结果,并生成VCF/BCF格式文件,BCF是VCF的二进制文件。再使用Bcftools对VCF/BCF格式文件进行SNP/Indel calling。其中,Bcftools是附属于samtools的程序。

    后面接Bcftools进行variation calling

    bcftools进行SNP calling:  www.jianshu.com/p/af93ca26ba44

    time $bcf mpileup ${precallbam} --fasta-ref ${reference} | $bcf call -mv -o ./out/bcfcall.vcf && echo " ** bcf snp calling done **"

    在进行SNP calling 时,必须选择一种算法,有两种calling算法可供选择,分别对应-c和-m参数。-c参数对应consensus-caller算法,-m参数对应multiallelic-caller算法,后者更适合多种allel和罕见变异的calling。

    -v参数也是常用参数,作用是只输出变异位点的信息,如果一个位点不是snp/indel, 不会输出。


    3.bam文件格式

    samtools view -H SRR1770413.sorted.markdup.bam

    详细讲解bam文件:https://www.jianshu.com/p/364e640d3c9f



    4.bftools合并vcf文件

    师兄用bftools合并vcf文件,这里的vcf文件不是gatk得到的,没有combine过。

    https://www.jianshu.com/p/59728f20d38c


    5 vcf 文件格式

    https://www.jianshu.com/p/38f734ae47f5

    https://www.jianshu.com/p/957efb50108f


    刚刚看到一个实战:

    https://www.jianshu.com/p/6f3198b7a070

    发现gatk在github上有workflow: https://github.com/gatk-workflows

    其中vcf文件生成步骤:

    1.HaplotyperCaller      2.MergeVcfs     3.GenomicsDBImport     4.GenotypeGVCFs

    在github上选择workflow打开,查看wdl文件中的command了解流程

    相关文章

      网友评论

        本文标题:10.18 bcftools下载使用 samtools语句

        本文链接:https://www.haomeiwen.com/subject/rzkjzftx.html